Post-Newtonian gravitational radiation and equations of motion 
via direct integration of the relaxed Einstein equations. 

I. Foundations 



Michael E. Pati and Clifford M. Will 

McDonnell Center for the Space Sciences, Department of Physics, 
Washington University, St. Louis, Missouri 63130 
(February 7, 2008) 



Abstract 

We present a self-contained framework called Direct Integration of the Re- 
laxed Einstein Equations (DIRE) for calculating equations of motion and 
gravitational radiation emission for isolated gravitating systems based on the 
post-Newtonian approximation. We cast the Einstein equations into their 
"relaxed" form of a flat-spacetime wave equation together with a harmonic 
gauge condition, and solve the equations formally as a retarded integral over 
the past null cone of the field point (chosen to be within the near zone when 
calculating equations of motion, and in the far zone when calculating grav- 
itational radiation). The "inner" part of this integral (within a sphere of 
radius 1Z ~ one gravitational wavelength) is approximated in a slow-motion 
expansion using standard techniques; the "outer" part, extending over the 
radiation zone, is evaluated using a null integration variable. We show gen- 
erally and explicitly that all contributions to the inner integrals that depend 
on 1Z cancel corresponding terms from the outer integrals, and that the outer 
integrals converge at infinity, subject only to reasonable assumptions about 
the past behavior of the source. The method cures defects that plagued 
previous "brute-force" slow-motion approaches to motion and gravitational 
radiation for isolated systems. We detail the procedure for iterating the solu- 
tions in a weak-field, slow-motion approximation, and derive expressions for 
the near-zone field through 3.5 post-Newtonian order in terms of Poisson-like 
potentials. 
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I. INTRODUCTION 



The motion of multiple, isolated bodies under their mutual gravitational attraction and 
the resulting emission of gravitational radiation is a long-standing problem that dates back 
to the first years following the publication of general relativity (GR). It has at times been 
controversial (for a thorough review see 0). In 1916 Einstein calculated the gravitational 
radiation emitted by a laboratory-scale object using the linearized version of GR ||. Some 
of his assumptions were questionable and his answer for the energy flux was off by a factor 
of two (an error pointed out by Eddington ||). In 1916, de Sitter f|] derived N-body 
equations of motion in what later would be termed the post-Newtonian (PN) approximation. 
However, his equations contained an error that was discovered in the course of a disputed 
claim by Levi Civita H that the center of mass of a binary star system would suffer a "self- 
acceleration" . Eddington and Clark || corrected the error, and found no self-acceleration. 
Einstein, Infeld and Hoffman J7| attempted to demonstrate explicitly that the Einstein 
equations alone imply equations of motion, by matching solutions of the vacuum equations, 
expanded in a weak-field, slow motion approximation, to fields representing the near-zone 
fields of "point" masses, working to first PN order. The result was the well-known EIH 
N-body equations of motion. Other highlights in this early history of the problem of motion 
include the development of the post-Newtonian approximation for fluid sytems by Fock || 
and Chandrasekhar ||, its extension by Chandrasekhar and later workers to 2.5PN order 
rnyi!] , and the development of equations of motion for spinning bodies by Papapetrou and 



collaborators ||12|| . 

Gravitational theory presents one problem essentially identical to that of electromagnetic 
theory: how to mesh the natural solution of the field equations in the near zone where the 
bodies reside, which involves slow-motion expansions and instantaneous fields, with the 
solution in the far zone, which involves retarded fields. Such a meshing is needed if one is 
to calculate the effects of gravitational radiation reaction that result from the emission of 
enegy and angular momentum to infinity. One approach to resolving this problem was that 
of matched asymptotic expansions. Although well-rooted in applied mathematics, it was first 
expounded in 1970 as a powerful technique for electromagnetic and gravitational problems 
by Burke |13| . Another, related approach is the "post-Minkowskian" framework, elaborated 
and developed most fully by Blanchet and Damour and their collaborators flH)-|T9"|l. 

A second important problem of gravitation, which distinguishes it from electromag- 
netism, is the non-linearity of Einstein's equations. Gravitation itself acts as a source of 
gravitation. Consequently this source extends over all space, resulting in the possibility of 
divergent or ill-defined integrals. In many ways, this has been the most serious difficulty 
to overcome. Techniques for resolving it have ranged from sweeping the difficulties under 
the rug, to the sophisticated analytic regularization methods of the post-Minkowskian pro- 
gram. A central thrust of this paper is to present a straightforward method for resolving 
this difficulty. 

A third "problem" , which is less a problem for gravitation than it is for electromagnetism, 
is that of "point" sources. In electromagnetic theory, where there is a belief that fundamental 
charges like the electron are point-like, the singular nature of the fields at the source has 
led to problems of mass regularization, especially in deriving equations of electromagnetic 
radiation reaction; it also raises issues of the boundary between classical and quantum 
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electrodynamics. In gravitation theory, this is less of an issue of principle, because the 
primary interest is in the motion of and radiation by astrophysical systems, whose members 
are clearly not point masses. Instead, the use of "point", i.e. delta-function sources is 
meant as an efficient means of approximating the mass distribution of bodies that are nearly 
spherical and that are small compared to the typical separation between them, so that tidal 
effects, which depend on the finite size of the bodies, can be ignored. Here the issue is how 
to make use of a point mass approximation (which simplifies many calculations) in a way 
that captures all the physics without introducing spurious effects. 

A fourth problem is of a technical nature: in electromagnetic theory, radiation damping 
in the equations of motion occurs at order (v/c) 3 beyond the simple Coulomb forces between 
charges, and is relatively easy to compute in a systematic approximation method, modulo 
the other problems listed above. By contrast, gravitational radiation damping occurs at 
order (v/c) 5 beyond Newtonian gravity, and requires a higher order of approximation that 
captures all relevant contributions. Over the years, numerous inequivalent results have been 
quoted for the leading gravitational radiation reaction effects. One finds published papers 
in which the coefficient in the relevant formula has ranged from —21/16 to the correct 
coefficient of unity); a study by Walker and Will [20] showed that the divergent results were 
all the simple consequence of missing one or more terms that contribute to the final answer. 

These four "problems" were the origin of the so-called "quadrupole controversy" , which 
arose from a critique by Ehlers and colleagues [21] of the foundations of the quadrupole 
formula for the leading-order gravitational radiation energy flux and orbital damping. This 
critique had the beneficial effect of spurring new research on those foundations, including 
a study of the systematic structure of the approximation sequence of Einstein's equations 
in a slow-motion, weak-field approach; analysis of energy balance as an argument for con- 
necting the far-zone energy flux to the near-zone damping forces, and elaboration of the 
post-Minkowskian approach, among others (see [|lj for a review). The work inspired by the 
Ehlers critique served to confirm the quadrupole formula and to strengthen its foundations. 
The ultimate test, of course, came in 1979 with the announcement of the measurement of or- 
bital damping of the binary pulsar PSR 1913+16 in agreement with the quadrupole formula 
|22|| ; current results agree to better than 0.5 percent p3| . 

The problem of motion and radiation has received renewed interest since 1990, with the 
proposal for large-scale laser interferometric gravitational-wave observatories, such as the 
LIGO project in the US, and the realization that a leading candidate source of detectable 
waves would be the radiation-reaction driven inspiral of a binary system of compact objects 
(neutron stars or black holes) |24|. Furthermore, it was noted [^5| that the leading method 
for data analysis of signals from such systems, optimal matched filtering, would require 
theoretical template waveforms that are accurate (primarily in the evolution of the orbital 
frequency or phase) well beyond the leading-order prediction of the quadrupole formula, 
possibly as high as corrections of order (v/c) 6 . 

This presented a major theoretical challenge: to calculate the motion and radiation to 
very high PN order, a formidable algebraic task, while addressing each of the problems listed 
above sufficiently well to ensure that the results were physically meaningful. This challenge 
was taken up by three groups of workers. 

One group, headed by Blanchet, Damour and Iyer [|14|-|T9||, used the post-Minkowskian 
(PM) approach to derive the gravitational waveform, equations of motion and energy flux 
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explicitly to 2PN order (0(v/c) A ) and beyond. The idea is to solve the vacuum Einstein 
equations in the radiation zone in an expansion in powers of Newton's constant G, and 
to express the asymptotic solutions in terms of a set of formal, time-dependent, symmetric 
and trace-free (STF) multipole moments [ 26| . Then, in a near zone within one characteristic 
wavelength of the radiation, the equations including the material source are solved in a slow- 
motion approximation (expansion in powers of 1/c) that yields both equations of motion for 
the source bodies, as well as a set of STF source multipole moments expressed as integrals 
over the "effective" source, including both matter and gravitational field contributions. The 
solutions involving the two sets of moments are then matched in an intermediate overlap 
zone, resulting in a connection between the formal radiative moments and the source mo- 
ments. The matching also provides a natural way, using analytic continuation, to regularize 
integrals involving the non-compact contributions of gravitational stress-energy, that might 
otherwise be divergent. 

The second group of Will, Wiseman and Pati use the approach described in the present 
paper, Direct Integration of the Relaxed Einstein Equations (DIRE), which builds upon 
earlier work by Epstein, Wagoner, Will and Wiseman pT|-|32|j. Like the PM approach, it 
involves rewriting the Einstein equations in their "relaxed" form, namely as an inhomoge- 
neous, flat-spacetime wave equation for a field h al3 , whose source consists of both the ma- 
terial stress-energy, and a "gravitational stress-energy" made up of all the terms non-linear 
in h a/3 . The wave equation is accompanied by a harmonic or deDonder gauge condition 
on h af3 , which serves to specify a coordinate system, and also imposes equations of motion 
on the sources. Unlike the post-Minkowskian approach, a single formal solution is written 
down, valid everywhere in spacetime. This formal solution, based on the flat-spacetime 
retarded Green function, is a retarded integral equation for h 01 ^ 3 , which is then iterated in 
a slow-motion (v/c < 1), weak-field (||^ Q/3 || < 1 ) approximation, that is very similar to 
the corresponding procedure in electromagnetism. However, because the integrand of this 
retarded integral is not compact by virtue of the non-linear field contributions, one quickly 
runs up against integrals that are not well defined, or worse, are divergent. Although at 
the lowest quadrupole and first PN order, various arguments were given to justify sweeping 
such problems under the rug p7|j28[] , they were not very rigorous, and provided no guaran- 
tee that the divergences would not become insurmountable at higher PN orders. Indeed it 
is straightforward to demonstrate that at second post-Newtonian (2PN) order, the rug is 
indeed pulled out from under such arguments. 

DIRE resolves these problems. The solution of the relaxed Einstein equation is a re- 
tarded integral, over the past null cone of the field point. The part of the integral that 
extends over the intersection between the past null cone and the material source and the 
near zone is approximated by a slow-motion expansion involving spatial integrals of mo- 
ments of the source, including the non-compact gravitational contributions, just as in the 
post-Minkowskian and Epstein- Wagoner frameworks. But instead of extending the spatial 
integrals to infinity as was implicit in earlier procedures, we terminate the integrals at the 
boundary of the near zone, chosen to be at a radius 71 given roughly by one wavelength 
of the gravitational radiation. For the integral over the rest of the past null cone exterior 
to the near zone ("radiation zone"), we use a change of integration variables to convert 
the integral into a convenient, easy-to-calculate form, that is manifestly convergent, subject 
only to reasonable assumptions about the past behavior of the source, that fully accounts for 
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the retardation of the fields comprising the source stress-energy, and that does not involve 
an explicit slow-motion expansion. This transformation was suggested by our earlier work 
on a non-linear gravitational-wave phenomenon called the Christodoulou memory [RO] (it is 
also implicit in Appendix D of |T^|). Not only are all integrations now explicitly finite and 



convergent, we can show explicitly that all contributions from the near- zone spatial inte- 
grals that depend upon the radius TZ are actually cancelled by corresponding terms from the 
radiation- zone integrals, for all powers of TZ (including In TZ), and for any order in the PN 
expansion. Thus the procedure, as expected, has no dependence on the arbitrarily chosen 
boundary radius TZ of the near-zone, and provides a simple practical method for regularizing 
integrals over non-compact sources. 

The ultimate products of this work will consist of equations of motion, gravitational 
waveforms, and energy flux expressions, in reasonably ready-to-use forms. The equations of 
motion for a binary system will have the schematic form 

d 2 x/dt 2 = -(Gmx/r 3 )[l + 0(e) + 0(e 2 ) + 0(e 5/2 ) + 0(e 3 ) + 0(e 7/2 ) + ...], (1.1) 

where m is the total mass of the binary system, x = xx — x 2 is the separation vector and 
r = |x|. The expansion parameter e is related to the orbital variables by e ~ Gm/rc 2 ~ 
(v/c) 2 , where v is the relative velocity. The leading term is Newtonian gravity. The next 
term 0(e) is the first post-Newtonian correction, which gives rise to such effects as the 
advance of the periastron. The terms of 0(e 2 ) and 0(e 3 ) are non-dissipative 2PN and 
3PN corrections. The 0(e 5 / 2 ) and 0(e 7 / 2 ) terms are the leading 2.5PN and post-Newtonian 
corrected 3.5PN gravitational radiation-reaction terms. (We do not include in this discussion 
contributions from spin, whose ordering in the PN hierarchy for compact bodies follows a 
special convention.) Explicit formulae for terms through various orders have been calculated 
by various authors: non-radiative terms through 2PN order [|i|,^3]-|36|| , radiation reaction 



terms at 2.5PN and 3.5PN order [37-39], and non-radiatve 3PN terms [4C-43 

In order to derive equations of motion to the 3.5PN order shown, one must derive the 
near-zone metric g a p as a function of spacetime and a functional of the source variables 
to 3.5PN order, which implies the following specific PN orders: goo through 0(e 9 / 2 ), goi 
through 0(e 4 ), gij through 0(e 7//2 ). In this paper we provide the required expressions in the 
form of (a) Poisson-like integrals of source densities, /(£, x')|x — x.'\ p d 3 x', where f(t, x') 
could be proportional to source stress-energy densities, and thus have compact support, or 
could be a function of other potentials, and thus extend over the entire near-zone region of 
integration A4; and (b) expressions involving time derivatives of source multipole moments 
Mijk... contracted with spatial vectors x l x^x .... These expressions can be simplified, iter- 
ated, and evaluated more explicitly, depending on the application envisioned ( "point" mass 
binary system, spinning masses, perfect fluid distributions, etc.). 

The second product will be expressions for the gravitational waveform, given schemati- 
cally by 

where [i is the reduced mass, and the subscript TT denotes the "transverse-traceless" part. 
The leading contribution Gfiv 2 /Rc 4: ~ GI^ /Rc A is the standard quadrupole formula. Ex- 
plicit formulae for all terms through 2.5PN order have been derived by various authors 



h ij = ^4 Pi 1 + °( eV2 ) + + °( e3/2 ) + + °^' 2 ) + ■ ■ \} TT , (1-2) 



29U32, 44-48 
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From the waveform, one can also derive expressions for fluxes of energy, angular momen- 
tum and linear momentum; the energy flux can be written in the schematic form 

dE/dt = (dE/dt) Q [l + 0(e) + 0(e 3/2 ) + 0(e 2 ) + 0(e 5/2 ) + 0(e 3 ) + ...], (1.3) 

where (dE/dt)g denotes the lowest-order quadrupole contribution, 

A third approach focusses on the limit in which one body is much less massive than the 
other, and employs black-hole perturbation theory to derive the gravitational waveform and 
energy flux, for particles orbiting both rotating and non-rotating holes. This method yields 
both numerically accurate results as well as analytic PN expansions up to orders as high 
as (f/c) 11 [^ , ^9H 53|. Work is currently in progress to extend these methods beyond the 



test-mass approximation, in an effort to compute corrections to first order in fi/M, the ratio 



of the mass of the particle to that of the black hole [p^-[56| . 

This is the first in a series of papers that will treat the problem of motion and gravi- 
tational radiation systematically using the DIRE approach. This paper lays out the foun- 
dations of the method, and derives formal solutions to the near-zone fields through 3.5PN 
order (order {v/c) 7 beyond Newtonian gravity), in a form useful for future applications. 
Subsequent papers in the series will derive the explicit equations of motion and near-zone 
gravitational fields for binary systems of compact object through 2PN order, and deal with 
radiation reaction at 2.5PN and 3.5PN order. 



Our conventions and notation generally follow those of [ 57| , |2"u ]. Henceforth we use units in 
which G = c = 1. Greek indices run over four spacetime values 0, 1, 2, 3, while Latin indices 
run over three spatial values 1, 2, 3; commas denote partial derivatives with respect to a 
chosen coordinate system, while semicolons denote covariant derivatives; repeated indices are 
summed over; rf v = r]^ u = diag(— 1, 1, 1, 1); g = det(g fMU ); = (a tj + a jt )/2; = (a lj — 
a- 7 *)/2; e 1 ^ is the totally antisymmetric Levi-Civita symbol (e 123 = +1). We use a multi-index 
notation for products of vector components and partial derivatives, and for multiple spatial 
indices: = x l x^ . ..x k , dij.„k = didj...dk, with a capital letter superscript denoting 

an abstract product of that dimensionality: x Q = x 11 x %2 . . .x lq and 8q = d^d^.-.d^. Also, 
for a tensor of rank Q, f® = f lll2 -- l i. Angular brackets around indices denote symmetric, 
trace-free (STF) combinations (see Appendix [A] for definitions). Spatial indices are freely 
raised and lowered with and 5^. 



II. FOUNDATIONS OF DIRE 
A. The relaxed Einstein equations 

We begin by reviewing the method for recasting the Einstein Equations 

R ap _ g <xP R =%nT a P , (2.1) 

into their "relaxed" form. Here R a/S and R are the Ricci tensor and scalar, respectively, g a/3 
is the spacetime metric and T a/3 is the stress-energy tensor of the matter. We define the 
potential 

h«P = r) a P - (-g) l/2 g aP , (2.2) 
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(see e.g. |[26|| ) and choose a particular coordinate system defined by the deDonder or harmonic 
gauge condition 

h^, =O. (2.3) 

With these definitions the Einstein equations (|2. 1|) take the form 

nh af> = -167ir a/3 , (2.4) 

where □ = —d 2 /dt 2 + V 2 is the flat-spacetime wave operator. The source on the right-hand 
side is given by the "effective" stress-energy pseudotensor 

T a P = (-g)T al3 + (levr)-^ , (2.5) 

where A af3 is the non-linear "field" contribution given by 

A<# = ltor(-g)tf L + {h^, v h^ >(t -h^,^ hn , (2.6) 

and t°^ L is the "Landau-Lifshitz" pseudotensor, given by 

lQn(-g)tf L = g x »g» p h a \y^ p + \g x ^h x \ p h?» , v - 2g pu g x ^\ p h^, x 

+ h 2 g aX g^ - g aP g^){2g vp g aT - g pa g UT )h^ . (2.7) 

o 

By virtue of the gauge condition (|2.3|), this source term satisfies the conservation law 

r a % = , (2.8) 

which is equivalent to the equation of motion of the matter 

T aP , p = 0. (2.9) 

Equation ( |2.4| ) is exact, and relies only on the assumption that spacetime can be covered 
by harmonic coordinates. It is called "relaxed" because it can be solved formally as a func- 
tional of source variables without specifying the motion of the source. Then, the harmonic 
gauge condition, Eq. ( |2.3| ) or the equations of motion, Eq. ( |2.9| ) are imposed to determine 
the metric as a function of spacetime. 

Notice that the "source" in Eq. ( |2.4| ) contains a gravitational part that depends explicitly 
on h al3 , the very quantity for which we are trying to solve. Also, we can expect r Q/3 , which 
depends on the fields h a P, to have infinite spatial extent. Indeed the very outgoing radiation 
that we hope to calculate, will, at some level of approximation, serve as a contribution to 
the source, thus generating an additional component of the radiation. 

Another complication in Eq. ( |2.4|) is that the second derivative term h al3 , pi/ h^ v in the 
source really "belongs" on the left-hand side with the other second derivative terms in the 
wave operator. This term modifies the propagation characteristics of the field from the flat- 
spacetime characteristics represented by the d'Alembertian operator to those of the true 
null cones of the curved spacetime around the source, which deviate from the flat null cones 



7 



of the harmonic coordinates. Nevertheless, the DIRE technique automatically recovers the 
leading manifestations of this effect, commonly known as "tails". 

The material will be modeled as perfect fluid, having stress-energy tensor 

T afB = {p + p^u 13 + pg a(3 , (2. 10) 

where p and p are the locally measured energy density and pressure, respectively, and u a 
is the four-velocity of an element of fluid. Until we begin to apply our results to specific 
physical situations, such as compact binary stars, we will have no need to specialize T al3 
further. 



B. Near-zone and radiation-zone 

We consider the material source to consist of a bound system of characteristic size S, 
with a suitably defined center of mass chosen to be at the origin of coordinates, X = 0. The 
source zone then consists of the world tube T = {x a \r < S, — oo < t < oo}. Outside T, 

The fluid is assumed to move with characteristic velocity ti C 1. The characteristic 
reduced wavelength of gravitational radiation, A = X/2tt ~ S/v = TZ serves to define the 
boundary of the near zone, defined to be the world tube V = {x a \r < TZ, — oo < t < oo 
}. Within the near zone, the gravitational fields can be treated as almost instantaneous 
functions of the source variables, i.e. retardation can be ignored or treated as a small 
perturbation of instantaneous solutions. For physical situations of interest, up to the point 
where the post-Newtonian approximation breaks down, TZ ^> S. The region exterior to the 
near zone is the radiation zone, r > TZ. 



The formal "solution" to Eq. ( |2.4j ) with an outgoing wave boundary condition can be 



written down in terms of the retarded, flat-space Green function: 

fc „ = 4 r ^TO'- f|+ |x-x1) dV , (2.11) 



but is really just a conversion of the differential equation ( |2.4|) to an integral equation. It 
represents an integration of r Q/3 /|x — x'| over the past harmonic null cone C emanating from 
the field point (t,x) (see Figs. 1 and 2). This past null cone intersects the world tube T> 
enclosing the near zone at the three-dimensional hypersurface M . Thus the integral of Eq. 
( [2 consists of two pieces, an integration over the hypersurface M , and an integration over 
the rest of the past null cone C — N '. Each of these integrations will be treated differently. 
We will also treat differently the two cases in which (a) the field point is outside the near 
zone, and (b) the field point is within the near zone (Fig. [[]). The former case will be 
relevant for calculating the gravitational-wave signal, while the latter will be important for 
calculating field contributions to r af3 that must be integrated over the near zone, as well as 
for calculating fields that enter the equations of motion for the source. 



S 



FIGURES 



Field point 
(far zone) 



Field point 
(near zone) 




FIG. 1. Past harmonic null cone C of the field point intersects the near zone T> in the hyper- 
surface M. Top: field point in the far zone; Bottom: field point in the near zone. Inner integrals 
are over the hypersurface jV, and outer integrals are over the remainder C — N of the null cone. 



C. Radiation-zone field point, inner integration 

For a field point in the radiation zone, and integration over the near zone (inner integral), 
we first carry out the t' integration in Eq. fl2.11| ), to obtain 



Z#(f,x)=4 



AT 







x — 


x' 


,x') 




X 


-x' 





d 3 x' . 



(2.12) 



Within the near zone, the spatial integration variable x' satisfies |x'| < 1Z < r, where the 
distance to the field point r = |x|. Expanding the x'-dependence in both occurrences of 
|x — x'| in the integrand in powers of |x'|/r, it is straightforward to show that 



h <*P 



(t,x) 



g=0 



where 



M af5Q [ 



u 



[ T a P( U ,X.')x' Q d 3 x' 
J M 



(2.13) 



(2.14) 



In Eqs. fl2.13|) and (|2.14|) , the index Q is a multi- index, such that 8q = d^d^ . . . di q and the 
superscript Q in M a ^ denotes ii«2 • • - iq, with summation over repeated indices assumed. 
The integrations in Eq. ( [2.14 ) are now over the hypersurface Ai, which is the intersection 
of the near-zone world-tube with the constant-time hypersurface = u = t — r. Roughly 
speaking, each term in the Taylor series is smaller than its predecessor by a factor of order 
d<1, provided we restrict attention to slow-motion sources. 
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Note that the field and source variables appearing in the integrand r a/3 are evaluated at 
the single retarded time u; however, because the field contributions to r a/3 fall off as some 
power of r, one can expect to encounter integrals that depend on positive powers of the 
radius TZ of the boundary of integration, especially in some of the higher-order moments. If 
this boundary were to be formally taken to oo (as has been the conventional approach in the 
past), these integrals would diverge. Instead we shall demonstrate (Sec. |TT| and Appendix 
|B|) that such ^-dependent effects are precisely cancelled by contributions from the "outer" 
integral. 

For the gravitational-wave signal, we need only to focus on the spatial components of 
h a/3 , and on the leading component in 1/R, where R is the distance to the detector. Using 
the fact that Uj = —N*, where N = x/i? denotes the observation direction, we obtain 



m>*) = p E / ^'(«,x')(N.xr^' + 0(r 2 ) . (2.15) 

it m=0 m\ at" 1 J m 



D. Radiation-zone field point, outer integration 

By making a change of integration variable from (r' , 9' , 0') to (u' , 9' , 0'), where 

t-u' = r'+ |x-x'| , (2.16) 
we can write the integral over the rest of the past null cone C — M in the form 

hf_ M {t^)=A du'l \ {U /[; X) [r'(u',n'))W, (2.17) 
J-oo Jc-N t — u' — n' • x 



where, from Eq. ( 2.16 ) 

r'(u', Q') = [(t - u'f - r 2 }/[2(t -u'-ti- x)] . (2.18) 

This change of variables represents an integration first over the two-dimensional intersection 
of the past null cone C with the future null cone t' = vl + r' emanating from the center of 
mass of the system at tcM = u' (Fig. 0), followed by the ^'-integration over all such future- 
directed cones, starting from the infinite past, and terminating in the cone emanating from 
the center of mass at time u, which is tangent to the past null cone of the observation point. 

For explicit calculations, it is useful to choose the field point x to be in the z-direction, 
so that n' • x = r cos9', and to write the outer integral in the form 

hf_Jt,x) = A du' d(j)' ' , [r'(u',n')] 2 dcos9' 

Ju-21Z JO Jl-a t — U' — n' • X 

ru-ITl r T al3 (n' 4- r' -X" /N | 

+4/ du' i [r'(u r , IT)] W , (2.19) 



t — v! — n' • x 



where 



a(u') = {u- u'){2r -2Tl + u- u')/2rK . (2.20) 
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The incomplete angular integration in the first integral of Eq. (|2.19|) reflects the fact that 
for u > u' > u — 21Z, the two-dimensional intersections meet the boundary of the near zone. 
For u' < u — 2TZ, the angular integration covers the full 4ir. Note that r a/3 contains only 
field contributions evaluated in the radiation zone; because they are themselves retarded, 
the "time dependence" u' + r' = t — |x — x'| is approximately constant over each angular 
integration, since it follows the hypersurface t — |x| = u = constant, and the dominant 
contribution to the fields comes from |x'| < 1Z. This allows a kind of slow-motion, multipole 
expansion to be exploited in evaluating these integrals, despite their range well outside the 
near zone. 




u'<u-2!^ u-2$^<u'<u u' = u 



FIG. 2. Change of variables for the outer integrals. Vertical line represents the material source 
world line. The variable u' is constant on the two-dimensional intersection between the past 
null cone of the field point and a future null cone from the center of mass of the system. Left: 
u' < u — 21Z, the two cones intersect fully outside the near zone, so the angular integrations are 
complete. Middle: u — 21Z < u' < u, angular integration terminates where the intersection between 
the two cones meets the boundary of the near zone. Right: v! = u, the upper limit of integration; 
the two cones are tangent to one another. 



E. Near-zone field point, inner integration 



In this case, in Eq. (|2.11|) , both x and x' are within the near zone, hence |x — x'| < 2TZ. 



Consequently, the variation in retarded time can be treated as a small perturbation, since 
r af3 varies on a time scale ~ 1Z. We therefore expand the retardation in powers of |x — x'|, 
to obtain 



/^(t,x)=4£-— / r^(t,x')|x-x'rW, (2.21) 

m =0 Ul JJVl 



where M. here denotes the intersection of the hypersurface t = constant with the near- 
zone world-tube. This version will be used for explicit calculations of the near-zone metric 
for use in the equations of motion. However an alternative formulation will be useful for 
studying the ^-dependence of the inner integrals; substituting the general Taylor expansion 
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x — x'| m 1 = £^L (— l) q (q\) 1 x < Q d > Q(r > m x ), where < (>) denotes the lesser (greater) of 
x| and |x'|, we obtain 

*x)=4E^E i # / /(t,^^ 1 )^. (2.22) 



F. Near-zone field point, outer integration 

The formulae from Section UP , such as ( 2.18 ) and (|2.20|) , carry over to this case with 
the result, 



u-2Tl+2r 



du' 



u-2Tl 



2tt 



1 T al3 (u' + r',x' 



+4 



-27?. 



r Q/3 (u' + r',x'; 
i — u' — n' • x 



l-a t — tt' — n' ■ X 

[r'{u' : Q!)\ 2 d 2 Q! . 



r'(u',Q')] 2 d cos 9' 



[2.23) 



Notice that the v! integration ends at u — 21Z + 2r rather than u because that corresponds 
to the last future null cone that intersects points in the far zone. 



G. Iteration of the relaxed Einstein equations 

Because the field h af3 appears in the source of the field equation, the usual method of 
solution is to iterate: substitute h al3 = in the right-hand side of Eq. ( |2.11| ) and solve 



for the first-iterated ih a/3 ; substitute that into Eq. ( [2.1 1|) and solve for the second- iterated 
2h a/3 , and so on (imposing the gauge condition Eq. ( |2~3l ) consistently at each order). The 
general sequence of iterations is shown shematically in Fig. |3|. The matter variables itia and 
the (N — l)-iterated field N-ih a/S are used to determine N-iT al3 and N-i-h alS . Eq. ( |2.11|) 



then yields ^h 013 as a function of spacetime and a functional of the matter variables. Then, 
if one wishes to determine the motion of the source, one substitutes ^h a ^ into the matter 
stress-energy tensor, and obtains the equations of motion from ArV i g(ArT a!/3 ) = where atV^ 
denotes the covariant derivative using the iVth iterated field. If one wishes to determine the 
iVth iterated gravitational field as a function of spacetime (i.e. with the matter variables 
determined as functions of spacetime to a consistent order), then one only needs to solve the 
equations of motion tv_iV ' ^{n-iT ) = 0, which are equivalent to the iV-th iterated gauge 
condition ^h a ^ \p = 0. 
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N _ iT aV = T aV (mA . Nlh aV ) 
^-iT a P = T a P( JV _ 1 r«P;^ 1 / ! «P) 



i 



v /t a P(m A ; x»)= I [^_iT aP ] Ret Ix-xT 1 ^ 





No 


N -»JV+1 





I 



N large enough? 



I 



Yes 



Far field or near field? 



Near 



Far 



r«P = r a P(m A ;^«P) 



I 



Get m^(xM-) from 



I 



Equation of motion 
of sources 



Get m^(xM-) from 



N-l 



I 



Compute gravitational 
waveform far from 



sources 



FIG. 3. Structure of iteration procedure 



H. General structure of the outer integrals 

At the first iteration, the solution is simply linearized general relativity. With q/i = 



substituted into the right- hand- side of Eq. (|2.11), the outer integrals vanish, and the inner 



integrals over the special relativistic T a/3 have compact support. There is no 7£-dependence 
in the integrals, trivially. For field points outside the source (|x| > |x'|), within both the 
near and far zones , the first-iterated \h a ^ takes the form of Eq. ( |2.13|) . Since M a @® is 
a function only of u = t — r, the spatial gradients 8q produce only unit radial vectors h\ 
powers of r, and retarded time derivatives of M a @Q . Products of n l can be grouped into 
symmetric trace-free (STF) products n <L> , which are analogous to Ylm (see Appendix |A| 
for useful formulae related to STF products). Thus, outside the source, ih alS can be written 
as a sequence of terms of the form 

x) = f B , L (u)n <L> r- B . (2.24) 

At the second iteration, in the far zone, T al3 = 0, and 1 A a ^(u' + r', x') consists of products 
of spatial and temporal derivatives of \h a ^{u' + r', x'). It therefore can also be expressed as 
a sequence of terms of the form 
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A a/ V + r', x') ~ f B A u >' <L>r '~ B ■ ( 2 - 25 ) 

Whenever the source at a given (N — 1) iteration takes this form, it is straightforward 
to evaluate the general form of the outer integrals for the iVth iterate. Defining the new 
variables £ = (t — u')/r — 1 + {u — u')/r, y — n • n' = cos 6, we find, from Eq. ( |2.18| ), 

r' = r(C 2 -l)/2(C-y). (2.26) 

Substituting Eqs. ( p.25|) and (|2.26|) into Eq. (|2.19|) , and changing to integration variables 
C, y and 4>, we obtain 

Nhf-M B ^ L = 2 (-) - <L> J_ x P Mdy j m ^ 2 _ V l )B „ 2 fBAu - r(C - 1)K , (2-27) 



where ((y) = z + ^z 2 — 2zy + 1, z = TZ/r, and Ph(y) is the Legendre polynomial. 

For far-zone field points, z < 1; Taylor expanding fs,L( u — r(£ — 1)) about w, we obtain, 
for 5 > 2, 



= f^^^En,,^)^^^, (2.28) 



g=0 

where the coefficients T> q B L (z) are given by 



n,L^ = ^ (C V 2 _ ^ AbAC, a)d( - E fc B,L (1 + 2z r—pT > (2-29) 
where 

A B)i (C,a) = |£_ a PL(y)(C- V)* -8 ^, (2.30a) 

«= (C-l)(C + l-2z)/2z, (2.30b) 

dk^l(C)/dC= k<&\), m>l; (2.30c) 

fcg> L (C) = A B>L (C 2)/(C 2 - I) 5 " 2 ■ (2.30d) 

The case B = 2 is special, and leads to the result 

N hf_^ L = ^ £ f 2 , L (u - s)Q L (l + ^ds + n< L > ± V\ L {zV d ^^ , (2.31) 



q=0 



where 



^w = 4r/ K "W-, k^*- (2 - 32) 

where Ql(v) is the Legendre function. 

Notice that, for B ^ 2, the outer integral returns a result of the same generic form as 
the input function. The case B = 2 returns terms with a logarithmic dependence on r (via 
the Ql's); terms of this form are called "tails". 
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Similarly, for field points in the near zone, z > 1, we Taylor expand Jb,l{u — r(( — 1)) 
about u + r = t, and obtain, for B > 2, 

where the coefficients £% l{ z ) are given by 

(-)Q r 2z+l C q 9 i , -n ( — l—2z) p 

^» = V-Li ( C2 -l)B- 2 ^(C^K-Efcfc P+1) (l + 2^) 1 p, j , (2.34) 
and, for B = 2, 

N hf_ N2 L = ^ / 2 , L ( U - S )Q L (l + £) rfs + £ ^.(^^(fi , (2.35) 
where 

f_)9+l f /-2z-l 1 /-2Z+1 rl-a P T (y) ) 

^w = L irU WMuLwL J0y) dv \- (2 36) 

Notice that, for near-zone field points, the functions /g,i are evaluated at the local time t, 
not retarded time u. 



I. Cancellation of TZ dependence 

It is evident that the inner integrals and outer integrals for the field h a ^ will separately 
depend upon the radius TZ of the boundary between the near zone and the far zone. But 
since each integral was simply a rewriting of a piece of the original integral, Eq. ( [2.11 ), 



which had no TZ dependence, it is equally evident that the separate ^-dependences must 



cancel between the inner and outer integrals. In [32|, referred to hereafter as WW, we 
demonstrated such a cancellation explicitly for contributions to the gravitational waveform 
at 2PN order that depended on positive powers of TZ. Here we demonstrate the cancellation 
generally, for both near-zone and far-zone field points, for arbitrary powers of TZ (including 
In TZ) and to an order of iteration sufficient for our purposes. 

The proof proceeds by induction. First, as we pointed out above, the first-iterated field 
ih a P is trivially independent of TZ. 

Secondly, we assume that the (N — l)-iterated field does not depend on TZ, i.e. that all 
^-dependence cancels at this order of iteration. We wish to demonstrate that this implies 
cancellation of ^-dependence in the iV-iterated field. The proof consists of considering the 
limiting behavior of the inner and outer integrals for N h a/3 in the vicinity of |x'| ^ TZ. Here 
T a/3 vanishes, and we only need to consider n-iA- , which is a functional of N-ih al3 . We 
have already seen that, in the far zone, N_iA al3 can be decomposed into terms of the form 
fB,L{u)h <L> r~ B . (We consider tail contributions with lnr dependence separately.) Since 
the N — 1 iterated field does not depend on TZ by assumption, continuity of the fields means 
that N_iA a/3 will have this same form just inside the near zone. Thus we will calculate 
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the limiting behavior of the inner integral of a term of this form as the integration variable 
r' — > 1Z from below, and compare its 1Z dependence with that of the outer integral of the 
same term. 

For far- zone field points, we must calculate the TZ dependence of the moments M a ^ , 
and substitute into Eq. ( |2.13| ); after considerable algebra (see Appendix |FJ), we obtain, for 



the limiting behavior of jy/j^f as the integration variable approaches TZ from inside, 



sh# BtL - (I) ' n< L > t^%zV d ^j^ , (2.37) 



where 

<? J max 



}m 2 2+j-B \i( q _ m _j + 2L)}\ (2L - j)\ 

u _ o m!( g -m-j + 2L + l)! [±( q - m - j\(L - j)\ 



f z 3+L-B +q - 1/{3 + L _ B + q _ j) 3 + L _ B + q _ j ^Q 

[lull 3 + L- B + q-j = . ^ ' ' 



where j max = lesser of {q — m,L}, and q — m — j = even integer > 0. Eq. ( p.37| ) is of 



the same form as the outer integral for far-zone field points, Eq. fl2.28|) . The coefficients 
T>% L (z) from the outer integrals are most easily evaluated using computer algebra methods 
(we calculated the coefficients using independent Maple and Mathematica programs); the 
result is, for each B, L and q, 

v b%(z) + i z -dependent part of V q BL (z)} = 0. (2.39) 

Thus the TZ dependence cancels term by term. Similar cancellation occurs for near-zone 
field points, as well as for the case where the integrand has r~ B In r dependence. Details are 
given in Appendix |B|. 

This cancellation, while inevitable, has practical consequences, in the following sense. In 
calculating the inner contributions to the fields, we must integrate over a finite hypersurface, 
Ai, sources that extend throughout M.. Consequently, any such integral will have terms 
that are independent of 71, as well as terms that depend on 1Z q or In 1Z. Because we know 
that all terms of the latter form cancel with contributions from the outer integrals in the 
final expression for the field, we can drop them in any individual result. Similarly, we can 
drop all ^-dependent terms that arise in any individual outer integral. This provides a kind 
of regularization of integrals, that cures the problem of divergent integrals that haunted 
earlier slow-motion methods. In fact, one can show that there is a close connection between 
this method of regularization and the method of analytic continuation used by Blanchet 



Thus, our procedure for determining the field is to determine separately h°ff and /i^-Af 
to a given PN order, keeping only ^-independent terms in each expression, then sum them 
to obtain 

H* = h% + hf_ M . (2.40) 
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III. WEAK FIELD, SLOW-MOTION APPROXIMATION 



We now turn to a discussion of the number N of iterations needed to derive equations 
of motion or gravitational waveforms of a desired accuracy, in a weak-field, slow-motion 
approximation. 

We assume that, for the fluid source, 

v 2 ~ m/S ~ p/p ~ e <C 1 , (3.1) 

where e will be used as an expansion parameter. But from the nature of the iteration 
procedure, it is evident that each iteration of the field introduces corrections of order m/S. 
In terms of e, m and S the equations of motion ( |1.1| ) can be rewritten schematically as 

dv/dt ~ (m/S 2 )[l + 0(e) + 0(e 2 ) + 0(e 5/2 ) + 0(e 3 ) + 0(e 7/2 ) + ...], (3.2) 

where the terms inside the square brackets represent the Newtonian, post-Newtonian, 2PN, 
2.5PN (radiation- reaction), 3PN, and 3.5PN (radiation-reaction) terms respectively. For a 
term of order e , the largest number of powers of m/S that can appear in it (including 
one power from the m/S 2 prefactor) is N + 1. The radiation reaction terms of order e N+l l 2 
must contain an odd number of velocities (in order to be odd under time-reversal), thus the 
maximum number of powers of m/S for them is also N + 1. Since one iteration gives the 
Newtonian potential, which yields the Newtonian equations of motion (N = 1), then, to 
obtain the 1PN terms (N = 2), one must have the second iterated field, to obtain the 2PN 
and 2.5PN terms (N = 3), one must have the third iterated field, while to obtain the 3PN 
and 3.5PN terms (N = 4), one must have the fourth iterated field. 

Similarly, to obtain a result for the waveform accurate to the order of the quadrupole 
formula, h ~ X 2 - 7 / R ~ (m/R)(v 2 + m/S) ~ e 2 (N = 2), the second-iterated field is needed. 
Note that the term m/S in X y arises through the use of the Newtonian equation of motion. 
Then, to obtain the 1PN, 2PN and 3PN corrections to the quadrupole approximation, the 
third, fourth, and fifth-iterated fields are needed, respectively. This would be an impossible 
task, if it weren't for the judicious use of the conservation law, Eq. ( |2.8| ). Consider for 
example, the source n-i t ^ of the iVth iterated gravitational- wave field at/^, Eq. ( |2.15 ), 



specifically the leading, m = term. The conservation law, Eq. ( |2.8| ), converts n~i t%j into 
two time derivatives of n^it 00 x 1 x^ (modulo total divergences). Because of the slow-motion 
approximation, two time derivatives increase the order by e, and thus, to sufficient accuracy, 
only jv_2T 00 is needed in practice. An important caveat to this is that the surface terms that 
arise from the total divergences and the outer integrals must formally be evaluated using 
the N — 1 expressions. However, in practice, these terms contribute at sufficiently high 
order that they can be treated without resort to explicit N — 1 expressions. Effectively, the 
burden of accuracy has been shifted from the iVth-iteration of the field, to the iV — 1-iterated 
equations of motion, which enter via the two time-derivatives, and which are needed anyway 
to evaluate the field as a function of spacetime. Thus, for iV = 2, the leading quadrupole 
approximation, only t 00 = p is needed, together with the Newtonian equations of motion. 
This circumstance is responsible for the prevalent, but erroneous view that linearized gravity 
(one iteration) suffices to derive the quadrupole formula. The formula so derived turns out 
to be "correct", but its foundation is not (see |58| for discussion). 
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Thus, in WW, to evaluate the 2PN waveforms (fourth iteration), only second- iterated 
fields were needed in the source terms. For 3PN waveforms, only third-iterated fields will 
be needed. 



IV. FORMAL STRUCTURE OF NEAR-ZONE FIELDS 

A. Metric and stress-energy pseudotensor in terms of the fields 

We begin by defining a simplified notation for the field h a ^: 

N = h 00 ^ 0(e) , 
K* = h 0i ~ 0(e 3/2 ) , 
B lj = h ij ~ 0(e 2 ) , 
B = h* = ^h" ~ 0(e 2 ) , (4.1) 

i 

where we show the leading order dependence on e in the near zone. To obtain the equations 
of motion to 3.5PN order, we need to determine the components of the physical metric to 
the following orders: g o to 0(e 9 / 2 ), g 0i to 0(e 4 ) , and g^j to 0(e 7//2 ). From the definition 
fl2.2Q , one can invert to find g a p in terms of h af3 . Expanding to the required order, we find, 

g 00 = -(1 --N+ -N 2 - -N 3 + — iV 4 ) + -B(l --N+ -N 2 ) 
yo ° v 2 8 16 128 ' 2 1 2 8 ' 

+I(B <J 'B ij " - ^B 2 ) + ^K j K j - ^NK j K j + 0(e 5 ) , (4.2a) 

g 0i = -K\l - l -N- l -B+ ^-N 2 ) - K j B ij + 0(e 9 / 2 ) , (4.2b) 

9ii = S ij (l + \n- l -N 2 + ^N 3 - \nB + l -K k K k ) 

+B ij - ^B5 ij - K l K j + -NB ij + 0(e 4 ) , (4.2c) 
(- g ) = 1 + N -B-NB + K { K { + 0(e 4 ) . (4.2d) 

Notice that, in order to find the metric g a p to the desired order, we must obtain N to 0(e 9 / 2 ), 
K l to 0(e 4 ), B lj to 0(e 7 / 2 ) and B to 0(e 9//2 ). In fact, because B contributes linearly to g 00 , 
we will treat B differently from B lJ . 

Using Eq. ( |4.2|) , we can express the matter stress-energy tensor T alS , Eq. ( |2.10| ), as a PN 



expansion. However, the details of such an expansion will depend on the basic variables used 
to characterize the matter. For example, to discuss the structure of a star in a PN expansion, 
it is convenient to use the mass-energy density p and pressure p, together with an equation of 
state. However, to discuss the motion of compact bodies in an effective "point-mass" limit, 
it is more convenient to ignore the pressure totally, and to use the so-called "conserved", 
or baryon density, p* = p^/—gu°. For now, we follow the convention of Damour et ai, and 
define the quantities 



a 


— rpOO _|_ rp 


a 1 


rpOi 




— T^i 
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1 



(4.3) 
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We will express various potentials formally in terms of these densities, and later make a PN 
expansion of them in terms of the densities most appropriate to the application. 

Substituting the formulae for h a/3 and g a p into Eqs. Q2.6| ) and ( p.7| ) for A alS , we obtain, 
to the required order, 

A oo = _I( V iV) 2 



+ \-N 2 -NN- 2N' k K k + -K id (ZK^ + K iJ ) 
18 2 

1„„ 7 



+k j Nj - B ij N> ij + -ViV ■ VB + -N(VN) 2 

4 8 

+ \x k ' j B jk + -B jk ' l (B jk ' L - 2B kl ' j ) + -NB - -(VB) 2 + -NN' j K j 
I 4 4 8 4 

+ l N d N ,k B jk _ l K i N ,k( 3K j,k + AK kj) _ I^ 2 (VA^) 2 1 + 0(pe 4 ) , (4.4a) 
8 2 8 J 



+ (jSfK* — NIC* — 2K k K i,k — B lm K i,lm _|_ j£k,l^jjU,k _|_ j^ik,l _ j^kl,i^ _|_ jy,k j^ik 

—-NB ,j ~ - -N' { B - NN' k (K k >* - K' l > k ) - -NNN' { + -KHVN) 2 - -K k N' k NA 
4 4 v ! A 8 V ; 4 J 

+0{pt 7 / 2 ) , (4.4b) 
A y = h^N^ - ^'(ViV) 2 } 

+ |2ir fe ' (i ir j) ' fc - K k,i K k ' j - K i>k K j ' k + 2N' {i K j) + -N' {i B> j) 

--N(N'' l N' 3 - -5 ll (VN) 2 ) - 5 ij (K l > k K [k $ + Ar> fc ;f fc + -N 2 + -ViV • VB)\ 
2 2 8 4 J 

+ \2K' l k j + B k(i (K j) > k - K Kj) ) - 2B ihk K k - NB ij - B ij ' lm B lm 

+ B ik ' l (B^ l ' k + B^ k ' 1 ) — 2B kl '( i B^ k ' 1 + - B kl,i B kl ^ — - B^ B^ 

2 4 

-N(2K k ' (i K j) ' k - K k,i K k,i - K i ' k K j > k ) + K k K k ' {i N' j) - 2NN' {i K j) - -NN' (i K j) 

--N' k N' {i B j)k - -NN' (i B' j) + -{VN) 2 B lj + -N 2 {N>' l N> j - -5 lj {VN) 2 ) 
2 2 8 4 2 

+ -5 ij [(VB) 2 + 2NB + 8K k ' l B kl + 4B kl ' m B km ' 1 - 2B kl ' m B kl ' m + 3NN 2 + 8NN> k K k 
8 

+%NK l ' k K [k ' l] - AK k N' l K Kl + 2NK k N' k + N' k N' l B kl + 2NVN ■ VB]} 
+0(pe 4 ) , (4.4c) 
A* = --(ViV) 2 



+ i K i,k K [k,i] _ N ,k^k _ _ VN . VB _ _j^2 + _jv(ViV) 2 } 

+ ^2K k K k - 2B' k K k + 3B kl K k ' 1 -NB + ^NB + ^(VB) 2 - B M B lrn 
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Q 1 1 

f}kl,m r>kl,m , ^ r>M,m r>km,l _ nr jy-l,k jAk,l\ _ _ AT,l r^k r^k,l , MM> k i{ k 

4 2 2 

+iiVA^' fc is: fc - -N' k N ,l B kl + iiVViV • VB + -(ViV) 2 ^ + -NA> 2 - -iV 2 (ViV) 2 l 
4 8 4 8 8 8 J 

+0(pe 4 ) , (4.4d) 

where an overdot denotes d/dt. In the above expressions, terms grouped within braces make 
leading contributions of the same order. For example, in A 00 , the three groupings correspond 
to 0(pe), 0(pe 2 ), and 0(pe 3 ), respectively. 



B. Source moments and other integral quantities 

Throughout our calculations, a number of integrals of the source stress-energy pseu- 
dotensor occur, for example, in the multipole expansions of Eq. ( 12.141) . It is useful to define 



and collect these quantities, and to discuss their properties. All integrals are carried out 
over a constant time (or constant retarded time) hypersurface Ai, within the near- zone. In 
general, these integrals will have TZ dependence, but, in line with the foregoing discussion, 
we shall consistently drop such terms. The relevant integrals are: 

P^ = M^°= [ T^°d 3 x, (4.5a) 

JM 

lQ = M 00Q = I r 00 x Q d 3 x, (4.5b) 

JM 

jiQ = e iab M 0baQ = £ ia6 j ^^fiQ^ ( 4 . 5c ) 

JM 

v ijabQ=f x l* T i]i3 x b]Q d 3 x _ ( 45d ) 
JM 

By making use of the equations of motion T af3 ^ = 0, we can transform some of these 
integrals into other forms, modulo surface integrals at the boundary dM. of the near zone. 
For example, 



- / T^ j d 2 S j 
JdM 

ji = I r 3h X a d 2 S j 

JdM 



T = P i - i r Qj x i d 2 S j . (4.6) 

JdM 

These identities express the conservation of total energy, momentum and angular momen- 
tum, and uniform center-of-mass motion, modulo a flux of gravitational radiation from the 
system. In calculations, the surface terms must be checked carefully to see if they make 
^-independent contributions to the order considered. For the most part, such surface terms 
turn out to make no contribution. 

Henceforth, we shall set X % = X 1 = 0, which amounts to attaching the origin of coordinates 
to the center of mass of the system. 

Other useful identities include: 
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M H = jy + I [ T lm ( x ij )j + f m0 x ij }d 2 S m , (4.7a) 

2 2 JdM 

1 2 

J^ijk _ _j-ijk _| e lk(ijl\j) 

6 3 

+- / \r lm (x ijk ) x + f m0 x ijk ]d 2 S m - - I [r l[k x i]] + r l[h x^\d 2 S l , (4.7b) 

6 JdM 3 JdM 

A/fijQ — \ J^jQ _i_ ^ mkx(i &m\j)k 2 ...k q ( ^ . r\\ , 8(g ~ jO T>U'( fc l fe 2-fcg) 

"(g+l)(g + 2) 1 + (g + 2) 6 J (symk.g j+ + 



(g + l)(g + 2) Jax 

2 : / [ r '[*i x *1i*2...fc 9 + 7 J[fci a j]<fc 2 ...fc,l d 2 5 ( k . Q) ( 4 _ 7c ) 

(g + 2) JaM ' 

M 0jQ = ^_jiQ Q_ e rnj(k 1 j m \k 2 ...k q ) + _J_ I T 0rn x jQ d 2 Sm fajft 

q+l q+l q + l JdM 

where the notation (sym k:Q) means symmetrize on the indices k\ through k q , and the 
superscript notation m l*-) means that only the indices following the vertical line are involved 
in symmetrization. 



C. Near-zone field expanded to 3.5 PN order 

We now carry out the explicit expansion of the near-zone field through 3.5PN order, 
beginning with the inner integral, Eq. ( |2.21| ), applying the above identities where possible. 
Inserting powers of e to indicate the leading order of each term, we obtain the result 

jV> = 4e / - 1 ' d 3 x' + 2e 2 d 2 / ^(t.x'JIx - x'ldV - -e 5 / 2 I kk (t) 
Jm |x — x'| J m 3 

+-e 3 <9 4 / t^xOIx-x'IW 
6 Jm 

1 f (5) (5) (5) 

_ A e 7/2 J ( 4x H + 2r 2 S kl^ jkl^ _ Ax k jkll^ + Jkkll^ 

30 

+ J-e 4 a 6 / r^^xOlx-x'lW 
180 t Jm v 71 1 

1 f (7) (7) (7) 

e 9/2 ^ 3r 4 J fcfc (t) +12rV j VHt) -12rV T kk (t) 

1260 1 w w w 

(7) (7) (7) (7) (7) 

+N 9M + 0(e 5 ) , (4.8a) 

i^ = 4e 3 / 2 / ^ (t,x/ W + 2e 5/2 <9i 2 / r°*(t, x')|x - x'|dV 
y V.M |x - x'| Jm 

o f (4) (4) ( 3 ) 1 1 f 

+ -e 3 ^ 3x fc X ifc (t) - l ikk {t) +2e mik J mk {t) \ + -e 7/2 <9 4 / r 0i (£,x')|x - x'| W 
9 6 Jm 
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1 



(6) (6) (6) (6) (6) 

^4 J on„2 J n-ik(±\ 1 n»,2 -rikkf±\ nn^,kl Tiki / +\ , i c~fe -rikll i ±\ o -rikkll i 



+— e 4 {30r¥r(t) -10r 2 r fefc (£) -20x M r fe '(t) +15x fc r fe "(t) -3J' ifcfe "(t) 



(5) (5) (5) (5) (5) 

20r 2 jT^t) +40x fc/ ,7 mfc (t) -lhx l J mkk {t) -30x fc J™* 1 ^) +12 J mlkk {t) 



+Kl )M + 0(e 9 / 2 ), (4.8b) 

fl y = 4e 2 /• ^yi d 3 x > _ 2e 5/2 jS (t) +2e 3gfl /" ^ ^ _ ^1^/ 

Jm |x — x'| Jm 

1 f (5) (5) (4) (3) 

--e 7/2 ^ 3r 2 r j '(t) -2x k l ijk {t) -8x k e mk{l J m]f} {t) +6 M l]kk {t) 

+ L* d f f r^'(t,x')|x-xfdV 
6 Jm 



1 



180 



(7) (7) (6) 

e 9/2 J 3r 4 jy( t ) _4rV J^(t) -16rVe mfc(i r'^t) 



(5) (5) (5) (5) 

+12r 2 M^' fcfc (t) +24x fc/ M ijk \t) -2Ax k M ijkl \t) +6 M ijkkll (t) 

+B% M +0(e 5 ). (4.8c) 

Explicit formulae for the boundary terms Nqm, K l dM and Bg M are given in Appendix 0. 
Through 3.5PN order, the terms in Eq. ( |4.8| ) divide naturally into two types: even terms, i.e. 
terms of integer powers in e in N and B^ and odd-half integer powers in K\ and odd terms, 
of odd-half integer powers in N and B^ and integer powers in K % . The even terms produce 
the leading Newtonian, PN, 2PN and 3PN contributions to the equations of motion, while 
the odd terms produce the gravitational radiation reaction forces. (Note that the even terms 
have odd contributions embedded within them, via contributions of the metric itself to r a/3 ). 
Through 3.5PN order, there is a clean division between even and odd terms, in the sense 
that even terms produce non-dissipative contributions to the equations of motion, while odd 
terms produce radiation reaction effects. At 4PN order this separation fails, because of the 
presence of tails - these are 0(e 3//2 ) modifications of the leading 2.5PN radiation-reaction 
terms, which result in disspative effects at 4PN order. We derive the leading contributions 
of these 4PN tail terms in Sec. |V1 C| . 

The outer integrals for near-zone field points turn out to contribute only beginning at 
3PN order (and, as we will see, do not contribute observable effects until 4PN order). This 
can be seen schematically as follows: for a source term of the form fB,L{u)h <L> r~ B , the 
outer integral has the form 

fit - r')(ti) <L> (r'r B r^ ~ r d -^r" l+l - B dT' ~ ^-M^ +2 ~ B , (4.9) 
M A 1 y ' |x-x'| Jn dtfl dti v ; 

where we have used the fact that |x| <C |x'|. The only possible ^-independent terms 
come from the case q = B — 2. Thus the outer integral gives a schematic contribution 
h>c-j\f ~ f^ B ~ 2 \t) where the superscript {B — 2) denotes B — 2 time derivatives. From Eq. 
( |4.4| a), the leading contribution to the source comes from (VA) 2 , where, from Eq. Q2.13Q , N 
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has the far-zone form N « 41/r + 2(3n <kl> l kl /r 3 + 3n <kl> l kl /r 2 + h kl l kl /r) + . . .. Taking the 
gradient of this expression and squaring, we get, schematically (ViV) 2 ~ X 2 /r 4 +I(I kl /r 6 + 
X kl /r 5 + 2 H /r 4 + . . .)• The first term (B = 4) gives no contribution, since X is constant 
to the order considered (X varies only via gravitational radiation energy loss). The second, 
third and fourth terms, (B = 6, 5, 4) together give h ~ II kk ^'{t). Since X kk ~ mr 2 , we 
find ft ~ (m/r) 2 i> 4 ~ 0(e 4 ), which is a 3PN contribution. Thus, for near zone field points, 
the outer integrals can be ignored until 3PN order. A similar argument for far-zone field 
points reveals that outer integrals begin to contribute only at 2PN order, as was found by 

ww m. 



D. Compendium of useful post-Newtonian near-zone potentials 



The even terms in Eq. (4Ji) have the form of ordinary Poisson-like potentials and their 
generalizations, sometimes called superpotentials. For a source /, we define the Poisson 
potential, superpotential, and superduperpotential to be 

P(J) = ^~f f^U V ' V 2 P(/) = -/, (4.10a) 

47T JM |X — X'| 

S{f) = ± [ _f(t,x')|x-x'|dV, V 2 S( /) = 2P(/) , (4.10b) 

47T JM 

SD{f) = ^f /(t,x')|x-x'|W, V 2 SD(f) = 12S(f). (4.10c) 

47T JM 

We also define potentials based on the "densities" a, a 1 and a 1 - 7 constructed from T a ^, 

£(/) = j M ^^^ d'x' = P(Anaf) , (4.11a) 

£*(/) - / t7t( |' x/)/( / |' X ° ^V = P(4™«/) , (4.iib) 

J.M x — x I 

e«(/)=/ a ^ x ' )/( f' x 'W = p ^ na i jf ) (411c) 

J.M |x — X'| 

along with the superpotentials 

X{f)= ( a(t,x.')f(t,x.')\x-x'\d 3 x' = S{4iraf), (4.12a) 

Y{f)= [ a{t,x')f{t,x')\x-x'\ 3 d 3 x' = SD(47Taf) , (4.12b) 

Z(/)=/ <x(i,x')/(£,x')|x-x'|W, (4.12c) 

and their obvious counterparts X*, X u , Y l , Y %3 , and so on. A number of potentials occur 
sufficiently frequently in the PN expansion that is it useful to define them specifically. First 
and foremost is the "Newtonian" potential, 

U= [ a{ t^\ d z x' = PUiza) = S(1). (4.13) 

JM X — X' 
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The potentials needed for the post-Newtonian limit are: 



$2 = 11(17), X = X(1). (4.14) 

Useful 2PN potentials include: 

= £*([/), =£(y*), 
X^X^l), X ij = X ij (l), 

X! = X«, X 2 = X(U), 

py = p(u>W) , p 2 = Pi? = $ 2 - ^t/ 2 , 

Gi = P(?7 2 ) , G 2 = P(C/C/) , 

G 3 = -P(U' k V k ) , G 4 = P{V^V^) , 

G 5 = -P(V k U' k ) , G Q = P(U' ij <5>\ j ) , 

G\ = P{U' k V k ^ + \P{U^U), H = P(U> ij P?) . (4.15) 

At 3PN order, the following potentials are useful: 

Y 1 = Y U , Y 2 = Y(U), Z = Z(1). (4.16) 

A variety of properties of these and general Poisson potentials are described in Appendix 
0. Note that, in evaluating Poisson potentials and superpotentials of sources that do not 
have compact support, our rule is to evaluate them on the finite, constant time hypersurface 
A4, and to discard all terms that depend on 1Z. 



V. EXPANSION OF NEAR-ZONE FIELDS TO 2.5PN ORDER 

We now turn to explicit evaluation of the near-zone fields and the metric to higher PN 
order, in terms of Poisson potentials and multipole moments. In addition to evaluating the 
inner integrals shown above, we must evaluate the outer integrals consistently at each PN 
order, to check whether any finite, ^-independent contributions result. 

In evaluating the contributions at each order, we shall use the following notation, 

N = e(N + eNi + e 3/2 A L5 + e 2 N 2 + e 5/2 A 2 . 5 + e 3 A 3 + e 7/2 A 3 . 5 ) + 0(e 5 ) , (5.1a) 

B = e\B x + e^B^ + eB 2 + e 3 / 2 77 2 . 5 + e 2 B 3 + e 5 / 2 B, 5 ) + 0(e 5 ) , (5.1c) 
fl« = e 2 (B l 2 J + e 1/2 B^ + eB l 3 j + e 3 / 2 ^) + 0(e 4 ) , (5. Id) 

where the subscript on each term indicates the level (1PN, 2PN, 2.5PN, etc.) of its leading 
contribution to the equations of motion. Notice that our separate treatment of B and B^ 
leads to the slightly awkward notational circumstance that, for example, B™ = B\. 
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A. The Newtonian and 1.5PN solution 



At lowest order in the PN expansion, we only need to evaluate r 00 = (— g)T 00 + O(pe) = 
<j + 0(pe) (recall that a 11 ~ eer). Since this has compact support, the outer integrals vanish, 
and we find 

N = 4U. (5.2) 

To this order, (-g) = 1 + 4U + 0{e 2 ). 

To the next PN order, we obtain, from Eqs. (|2.5|) , (|4.4j ) and ( |5.2|) , 

r 00 = a - a u + AaU - ^Vf/ 2 + 0{pe 2 ) , 

T°* = a* + 0(p e 3 / 2 ), 

r u = a u - — W 2 + 0{pe 2 ) , 
8n 

T v = O(pe) . (5.3) 

Substituting into Eqs. fl4.8|), and calculating terms through 1.5PN order (e.g. 0(e 5//2 ) in 
N), we obtain 

(5.4a) 
(5.4b) 
(5.4c) 

(5.4d) 
(5.4e) 

It is straightforward in this case to show that the outer integrals and surface terms give no 
^-independent terms. It is useful to illustrate the cancellation of an ^-dependent term in 
this simple case. In the far zone to Newtonian order, the field, from Eq. (|2.13| ), is given 
by N m 4X/r, where we focus on the monopole contribution. This contributes to A 00 in 
the far- zone a term of the form A 00 = — 14X 2 /r 4 . To evaluate the near- zone contribution 
of the outer integral of this term, we must evaluate the coefficient £% L in Eq. ( [2.33D with 
q = (no time derivative, since I is constant, to lowest order), B — 4, L — 0. From Eqs. 
( p.30| ) and (|2.34j ), this yields a contribution to N given by Nc~n = —7X 2 /TZ 2 . However, in 



iVi 


= 7U 2 - 4$! - 


f 2$ 2 + 2X , 


K\ 


= AV\ 




Bi 


= U 2 + 4$! - 


2$ 2 , 




o (3) 

= -- l kk (t) , 
3 W 




Bi.5 


(3) 

= -2 l kk (t) . 





evaluating N^, we encounter the Poisson potential — 14P(VU 2 ) = — 14P 2 (see Eq. ( |4.15| )). 
Upon integrating by parts and keeping the surface term at 1Z (see Eq. (|D3|a)), this gives a 
contribution 7U 2 — 14$ 2 + 7X 2 /TZ 2 , whose ^-dependent term cancels that from the outer 
integral. 

The physical metric to 1.5PN order is then 



4 



(3) 



goo = -1 + 2U - 2U 2 + X - ~ l kk {t) +0(e 3 ) , (5.5a) 

o 

^ = _4^ + 0(e 5 / 2 ), (5.5b) 

g i3 =5 lj {l + 2U)+0{e 2 ). (5.5c) 
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Notice that, in our formulation, the potential U is not a retarded potential; the retar- 
dation is expressed by the PN potential X and the 1.5PN term — |Z fcfc ( 3 )(t). This contrasts 
with the PM approach, where retarded, rather than Poisson potentials are used, and the 
retardation is expanded only much later in the computation. The apparently 1.5PN term 
—~Z kk ^(i) in goo actually doesn't contribute to the equations of motion at this order because 
it is purely a function of time, and the leading contribution is through a spatial gradient. 
As a consequence, the lowest-order observable contribution to radiation reaction is at 2.5PN 
order. (An alternative way to treat this 1.5PN term would be to absorb it in a redefinition 
of the time coordinate.) 



We note here the useful identity, which follows from Eqs. (|5.3| ), a = r + r 
(l/27r)V 2 (L r2 ) + 0(pe 2 ), whose consequence is, 



a(t, x)<f x = 1 + -V 1 + 0(le 2 ) , 
M 2 



where surface terms make no "^-independent contribution. 



B. The 2.5PN solution 



(5.6) 



At 2PN and 2.5PN order, we obtain from Eqs. (pip , |gb,c), ( f^ and flj[ 



r« = a ij + —{U' l U' j - -5 ij VU 2 ) + 0{pe 2 ) , 
Aix 2 

T 0i = a l + Aa l U + -U*VM + —UU'* + 0{pe 5/2 ) 
ix An 



(5.7a) 
(5.7b) 



Including outer integrals and boundary terms (which contribute nothing), we obtain, from 
Eq. 



B% = 4$1 J + AP l 2 j - 5 ij (2<S> 2 - U 2 ) , 

K\ = 8V* - 8¥ 2 + 8UV i + 16G" 7 + 2X l , 

(3) 

o ( 4 ) o (4) 4 {3) 

K l 2 5 = -x k l ik {t) — l ikk (t) + - e mik J mk (t) 
3 9 9 



(5.8a) 
(5.8b) 

(5.8c) 

(5.8d) 



The solutions for B 2 and B 25 , along with the earlier 1.5PN solutions must now be substi- 
tuted into (— g)T al3 and Eq. (|4.4|a,d), with the result 



r 00 = a - a 11 + AaU VU 2 

8n 

+a(7U 2 - 8$i + 2$ 2 + 2X) - 4a u U 

+— I^U 2 - AUU - 8U' k V k + 2V i ' j (3V j ' i + V id ) + 4V j U' j - 4U' ij ^ 

+8W ■ V$i - 4W ■ V$ 2 - -VU ■ VX - 10f/W 2 - 4C/«P 2 ij 

2 
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A (3) -i (3) 

+ aX kk (t)+—U^T^t), (5.9a) 

r" = <r u - — W 2 
8tt 

+4a ii E/ - — ( III 2 + AV id V [id] + W J U' j + IvU ■ Vxl . (5.9b) 
4vr 12 2 J v ' 

Substituting into Eq. (f4.8|a) and ( |4.8| c) and evaluating terms through 0(e 7//2 ), and verifying 
that the outer integrals and surface terms make no ^-independent contributions, we obtain, 

20 

N 2 = -16C/$! + 8C/$ 2 + 7UX + —U 3 - AV l V l - 16£($i) + E(X) + 8T, l (V l ) 

3 

1 (4) 

-2X 1 + X 2 + -Y -4Gi - 16G 2 + 32G 3 + 24G 4 - 16G 5 - 16G 6 - 16i? , (5.10a) 
6 

B 2 = UX + 4VV - S(X) - 8S l (^) + 16E"(C/) 

+2Xi - X 2 - 20Gi + 8G 4 + I6G5 , (5. 10b) 

1 (5) o (5) 1 (5) 

N 2 5 = (2x kl + r 2 5 kl ) l kl (t) +—x k l kll (t) -— I kkll (t) 

15 V ; w 15 W 30 W 

1 « (3) (3) 

+ ™U l kk {t) -AX' kl X k \t) , (5.10c) 
3 

1 (5) 9 (5) 8 (4) 2 (3) 

S 2 . 5 = — r 2 V\t) +-x k T ik (t) +-x k e mki J mi (t) — M ilkk (t) . (5.10d) 
3 9 9 3 



C. The far-zone field to 1.5PN order 

In anticipation of finding non-zero outer-integral contributions to the near-zone field at 
3PN order, we must determine the far-zone field to an order needed for the source A a/3 . Our 
foregoing discussion indicates that counting PN orders for outer integrals is different than 
the standard method, because the inverse radial variable r _1 < TZ^ 1 ~ v/S; in other words, 
when considering contributions to the outer integrals, additional powers of r in a term in 
the far-zone field can be regarded as increasing the effective order of that term by half a 
power in e. For example, expanding hffi = N_\f in the far-zone, Eq. ( |2.13| ), we obtain 



N N = A-+ l -d kl { ^) - ~d klm I )+■■■}, I-)-! I ) 




v r 2 \ r 
e e 2 

where the effective PN order of each term is indicated. In ordinary applications, the second 
potential in Eq. ( |5.11| ) would contribute a term of order e of the form n <kl> X kl /r 3 , which is 
simply the Newtonian quadrupole potential. But in the outer integral, this term contributes 
an ^-independent term only through several time derivatives, and thus its effective contri- 
bution is higher order, in fact of the same order as that of the term n kl I kl /r, which also 
comes from the second potential. 
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At this order, we must also be careful to include any outer integral and boundary con- 
tributions to the far-zone field. From the lowest-order far-zone field, we find, to the order 
needed, that A 00 = -14(X/r 2 ) 2 , A lj = A(X/r 2 ) 2 (n <ij> - 5 ij /6). Evaluating the coefficients 
V\ 2 and 2?2,o> Ec l- d!iD> we obtain, in the far zone, N C _ N = 7(J/r) 2 and B^_ N = (l/r) 2 h ij . 
Combining the multipole expansions of Eq. ( |2.13| ) with the outer integral contributions, we 
obtain in the far-zone, to the order needed, 



N = 4- + 2d kl - \d klm [l^m +7 % + 0(e 3 ) , (5.12a) 

r \ r I S \ r I r z 



K - = - 2dt (EM) + + l du + ^ (?M ) + 0(f 3) , (5 . 12b) 



B , = 2 EM + l Bk (E!M\ + l^ dk (^M\ + + 0(£ 3 } . (5 . 12c) 



It will turn out, however, that, despite the formal possibility of 3PN contributions from the 
outer integrals, the actual contributions will not begin until 4PN order (see Sec. [VI Cj) 



VI. EXPANSION OF NEAR-ZONE FIELDS TO 3.5PN ORDER 
A. BV and K* to 3PN and 3.5PN order 



At 3PN and 3.5PN order, we obtain from Eqs. Q, (fO]b,c), Q and (0) 



r ij = a 11 + — (U'W - -5 ij VU 2 ) 
4tt 2 ; 



+4<r <J 'Z7 + — \U' {l X' j) - lQV [i ' k] V lj ' k] + 8U' {l V j) 

47T L 

-6 ij (-X7U • VX - W [lM V m + 4U+ k V k) + -U 2 )) + 0(pe 3 ) , (6.1a) 
2 2 J 

r 0i = (x 4 + 4(j i C/ + + —UU' 1 

+a i (7U 2 - 8$i + 2$ 2 + 2X) 

+— \64U' k (vl k,i] - $!> M ) + 32UU' k V k ' i - WUU^V^ + 16U' i U' k V k 
IQtc l 

-2AV l (VU) 2 + lQU' k X [k ' l] + 128U' k Gf' l] - 32<f>fV [k ' i] - ld&JfV^ 

(3) 

-16X' k V [kA - 16t/$f + AWUU^ + Gf/X* + 6U' { X -16C/'^i + 16UV* 
-1GUV* - 32V k V i ' k - 16V*' w ($f + P 2 fe/ ) + 16£/' fc (<if + P l 2 k ) 
+16V k ' l ($f k + $f - $1^) + 16V k ' l (P?' k + P2' 1 - P^)} 

A (3) 1 f (3) (4) 

+-<t* I kk (t) +— V iM l k \t) -U> k l ik {t) , (6.1b) 
3 2n V / 

where the first line in each expression is the contribution through 2PN order obtained earlier. 
Substituting into Eq. (|4.8|b,c) and keeping contributions through 0(e 7 / 2 ), and checking that 
surface terms and outer integrals make no contribution to this order, we obtain, 
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Bl j = 16E' J '(P) + \ril \X j: ) - QAP(V [i ' k] V [j ' k] ) + 32P(U' (i V j) ) + 2X ij + 2S(U> i U> j ) 

+5 i3 (UX - AV k V k - E(X) + 8Z k (V k ) -X 2 - 8d - 8G 4 + 16G 5 ) , (6.2a) 

K\ = \2V' 1 V i + 16UV* - 16P$ 2 + AUP + 32UG\ + AV l X - 8$iV* + 8<& 2 V i 

_gy fe $f - ST^Pf - 16E(VJ) + 16E($ l 2 ) - 16E(?7F i ) - 4E(X*) - 32E(G* 7 ) 
-24E*($i) + 4E*(X) + 8E fcfc (V^) + 8E fe ($f ) + 8E fc (P 2 ifc ) + 8E ife (F fe ) 
+2AP(UUU' i ) + 2AP{U' k U' l V k ) + 32P([/'V 2 fc,i ) - 32P(P' fe $ 2 ' < ) + 64P(P' fc G 7 : ' i ) 

(3) 

+8P(P' fc X fc ' i ) + 16P(P' fe $f ) + 16P(P' fc Pf ) - 16P(f/' i <i> 1 ) + 6P(C/^ X) 
-16P(f/$0 + 6P(PX*) + 32P{UU' k V i ' k ) - 16P(UV*) + 16P(y fc,, $?' fc ) 
-16P(V rfc ' / $f *) + 16P{V k ' l P 2 ' k ) - 32P{V k V i ' k ) - 16P(V k ' l <S>' k ) + 8P(V k ' i X' k ) 

_i6P(y fc ' , p 2 w ' i ) - i6P(y i '' m $' 1 m ) - i6P(y i - ,m p.j m ) 

1 (4) 

+4X\U) - 4X(V™) + 8S(U> k V k ' 1 ) + 6S(P' 4 P) + - Y l , (6.2b) 

6 

i (5) o (5) 8 (4) 2 (3) 

= — r 2 P J (t) r j7c (t) +-x fe e mfc(i J m W(t) — M ijkk (t) , (6.2c) 

3 9 9 3 

1 (6) i (6) ? (6) -, (6) 1 (6) 

5 = _!_rV T k {t) --^r 2 T kk {t) -^x kl V k \t) + ^x k T kll {t) P fcfc "(t) 



' o ( 5 ) 4 ( 5 ) 1 ( 5 ) 1 ( 5 ) 9 (5) 

— r 2 ,7 m '(t) +— x M ,7 mfc (t) -— x l J mkk {t) -—x k J mk \t) +— J mlkk (t) 



l fi (3) (3) (4) 

+ —V i X kk {t) -AX iM Z k \t) +AX' k l ik {t) . (6.2d) 
3 



B. N and B to 3PN and 3.5PN order 

The expressions for r 00 and t u to 3PN and 3.5PN order are too lengthy to be reproduced 
explicitly. Instead, by substituting the expansions (|5.1|) into Eqs. ( |2.5|) and (|4.4| b,c), and 
keeping terms of 0(pe 3 ) and 0(pe 3,5 ), we obtain the formal contributions 

r 3 00 = a(N 2 -B 2 - N Q B 1 + K\K\) - a«(Ni - Pi) 

+T^(-^( 2ViV o • ViV 2 + (ViVx) 2 ) + ^-NqN, - iVoiV! - iVoiV! 
167T 18 4 

-2N$Kl - 2N(K{ + /\|-'(3/\.f + K^') + A^if* + A^iq - X<f P 3 j 

—N{ ij B 2 j + ](VN • VP 2 + VN X ■ VPi) + In^VNo) 2 
4 8 

+ 7 -N (VN ■ ViVO + A';-'/);' + \Bf k {Bf k - 23**) + -N^B X - ^(VP^ 2 

+ ^iV iV^i + - \n^K{{AK^ + 3#f ) - ^iV 2 (ViV ) 2 } , (6.3a) 

r 3 00 5 = v(N 2 . 5 - P 2 . 5 - N B 1£ ) - - P L5 ) 

+T^(-7 VAr o • ViV 2 . 5 + ^A> A>i. 5 - iVoJVi.5 - NoNx.5 - 2N^ S + N$Ki 5 



16vr I 4 4 
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+KWK& + KH) - N^Bf, - NfB'l + ±VN • YIL 

1 • 



+-iV 1 . 5 (ViV ) 2 + K[ J IV- + -NoBu + -.V U '.V,//^} , (6.3b) 
rf = ^(N, - B,) 

+T^(-i(2ViVo • ViV 2 + (ViVx) 2 ) + -IKI 'KJ 1 - N^K l 2 - N(k{ 

l07T 18 

-i(VJV • VP 2 + ViVx • VBi) - ^N Ni + ^N^VNq) 2 + ^N (VN ■ ViVi) 

+2iqiq - 2Si*iq + - N B\ + -iVo^i + -{VBrf - B[ ij B$ 

4 8 

+^Bf k Bf k + -BpB l 2 k ' j - .XjK^h'i'-' - ^N$K{Ki* + N N^k\ + ^N N$Ki 
-±NgNgB% + \n VN • VBi + ^(ViVo) 2 ^ + jjtf JV* - jjiV 2 (ViV ) 2 } , (6.3c) 

^5 = ^(^1.5-^1.5) 

+ ik{~l VN ° ' ViV2 - 5 + 2K ^ K ^ - N ^ - i ViV ° ' V£?2 - 5 

-|jVoiNTi. 6 + ^iV L5 (ViV ) 2 + :\B'1:Jk'; j - N B\. 5 

+In B 1j& - B?B% - \NiN$B% + ^(ViVo) 2 !^} . (6.3d) 

We have simplified the expressions slightly by taking into account the fact that A^ L5 , Bi 5 
and B % 2s are purely functions of time, so that spatial gradients of them vanish. To obtain the 
full expressions, one substitutes for N , N%, B 1 , K[, etc. from Eqs. (|5.2j ), (|5.4j) , (|5.8j) , ( p.lO| ) , 
and (|6.2j ). Substituting this into Eq. (|4.8[a) and (|4.8| c) (the latter contracted on indices ij), 
and including surface terms and outer integrals, we obtain the final 3PN and 3.5PN results 
for N and B: 
iq 

N 3 = — U A - 28f/ 2 $i + 14[/ 2 $ 2 + 10U 2 X - WX X + AUX 2 + AUE(X) 
6 

-AUG! - 56UG 2 + 112UG 3 + 80UG 4 - 6AUG 5 - 56UG 6 - 56UH + ^-U Y 
+32UE i (V i ) - 56C/E($i) - SUYf^U) + 10$ 2 - 8$i$ 2 + 2$ 2 - 8$iX + 4$ 2 X 
-16\rV$ + 16V% 2 - 32V l G\ - AV i X i + -X 2 - 2$\ j ^ - A^P? - 2P l 2 j P l 2 ] 
-8£(f/$i) + 36E(Gi) - 8E(G 2 ) + 16E(G 3 ) - 48E(G 4 ) - 8E(G 6 ) - 8E(P) 

1 (4) . . . . . . . . .... 

-8£(Xi) + y^S(F) + 16E 4 (V 2 4 ) - 16E l ($ l 2 ) - 16E l (£A^) + 32E i (G 7 ) + AY?(X l ) 

+8E l ' l '([/ 2 ) + ^E^O + 4£ ij ($i i ) + !>;''( /^) - 8E(E($i)) + 64E(E*(1/ 1 ')) 
-56E(E"([/)) - 32P{U 2 U) - 28P(UU 2 ) + 16P(tf$i) + 16P(?7$i) 

(4) . (3) 

+16P(£/$j) - 8P(£/$ 2 ) - 8P(t/ X) - 4P(£7 X) - 8P(UX) + 32P(UU' l V l ) 
-64P(f/f>' i 1/ i ) - AOPiUU^) + 32/ , (ri' ; '\'' r ) + 16P(£A^V J>i ) 

(3) 

+32P(U' l V 2 l ) - 32P(U' l & 2 ) + 64P(fW 7 ) + 8P(f/' 1 X J ) 
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+6AP(U'W 2 ) - 64P(f/'*V 2 i ) - 128P(U^G\) - 16P(C/'0P) 
+AP(U> i U> j <S>\ j ) + AP(U A U' J P l 2 j ) - QAP(U' ij i: ij (U)) - 8P(U' ij S^U' 4 )) 
-128P(U' ij P(U' {i V j) )) - lQP(U' ij P(U' {i X' j) )) + 25QP(U' ij P(V [k ' i] V [k ' j] )) 



(3) 



-8P(U' ij X ij ) + 32P{V i ^>t) - 16P(l /i $2 i ) - IQP{V 1 X'*) + \QP{V^ j ^) 
+16P(V i ' j P?) + 96P(\/ ij \/ 2 j ' i ) - 96P(V l ' j &/) + mPiV^G?) 
+2AP(V i ' j X j ' i ) + 8P{V i V i ) - 16P(\> i $0 + 8P{V i X> i ) - 8P($f V^) 
-16P($pPf '*) - 8P(P 2 ij ' k Pi h ' i ) + 16P($f + 16P($f P 2 ij ) - 8P(<S>\ j X' ij ) 

.... .... 1 (4) 1 (4) 1 (6) 

-8P(*W) - 8P($f P«) - 8P(P 2 ^) - - +- F 2 +— Z 

-8X($i) + ^X(X) + AX 1 ^) - 25'(f7 2 ) - 8S(C/C/) 
8S(U' ij $[ j ) - n.SS /•••'/';') + 12S(l^V^) + 8S(£W l ) - 16,S(\/ l f7' 1 ) , (6.4a) 
B 3 = \\J^ - 4f/ 2 $ 1 - 2f/ 2 $ 2 + WUV^ - \2UGx - 8UG 2 + 16UG 3 + I6C/G4 

1 (4) 

-8C/G 6 - 8C/# + —U Y -8C/E($i) + 8[/E M (f/) - 2$ 2 + 8$i$ 2 + 2$ 2 

+ 16 W 2 l - 16V*$* 2 + 32V l G\ + 4\/ l X* + ^X 2 - 6<S>\ j <S>¥ - \2^P l 2 j - 6P 2 iJ f P 2 ij 
+8E(C/$i) + 12E(G 1 ) + 8E(G 2 ) - 16E(G 3 ) - 16E(G 4 ) + 8E(G 6 ) + 812(H) 

1 (4) 

E(y) - 16E l (V2 l ) + 16E l ($ 2 ) - 48E*(fA/ i ) - 32E i (G 7 ) - 4E*(X l ) + 24E M (f/ 2 ) 

-28E"($i) + 8£"(X) + 12E^($i ) + 12E^(P 2 ij ) + 8E(E($i)) - 8E(E"(C/)) 

. (3) ... 

-68P(f/f/ 2 ) - 16P(C/$i) + 8P(C/$ 2 ) + 48P(f/$i) - 20P(C/ X) - 56P(C/C/' l y l ) 
-32P([/CW*) + 32P(E/V* J V* J ') - 16P(C/y ij V^) - 32P(U' i V 2 i ) + 32P([/' i $ 2 ) 



-64P(CW 7 ) - 8P(U> i X 1 ) - AP(U^U^^I) - AP{U' l U' j P' 2 j ) - 32P(V i <S>() 
+16P{V i <b 2 i ) + 48P(l^'$i j ) + A8P{V hj P 2 ] ) + 32P(V i ' j V 2 j ' i ) - 32P(V^'$j i ) 
+64P(\/ ij 'G , f) + 8P(V i ' j X j ' i ) + 2AP(V i V i ) + 16P(V^) - 8P{V i X> i ) 
+8P($f fc $f >*) + 16P($f fe Pf '*) + 8P(P 2 ij ' k Pi k '*) - 16P($f - 16P($f P 2 ij ) 



(3) 




(6.4b) 
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/16 34 16 \ ^ 1 ^ 

+ — r 2 U x l X' 1 + —X VHt) (IQx'U - 52X> 1 ) V 33 (t) 

V15 15 5 / 45 

/A 2 \ (5) 2 (5) 

+ -x k X' ij - —Y> ijk i i3k (t) -—u r ijj (t) 



9 27 / w 15 

4 . . 8 ■ • • 10 • A (4) • (4) 8 • ■ (4) .. 

- ( -X 1 ' 3 - -x l X> 3 + — Y^ T 3 (t) +8X V 3 (t) — X' 1 T 33 (t) 
v 3 3 9 / 9 

.. 2- \ (3) 
UUX> %3 + 2S(X«) - 4Xf + 2Xg 3 + 

1 (3) o (3) 4 (3) 

+-(70f/ 2 - 16$i + 20$ 2 + 4X) X"{t) +-U M 33kk (t) — X' ij M ijkk (t) 

3 3 3 

16 (4) 4 (4) 16 (3) 

+— x k X' i3 e mk{i J ml3 \t) — {8x k U - 5X' k )e mkj J m3 \t) + —X' k e mkj J mj {t) , (6.4c) 
9 9 9 

1 (7) i (7) 1 (5) 1 (5) 

B 3 5 = r 4 J^'(t) +— r 2 ^ r^'(t) — (Y' ij - 6x l X> 3 ) V 3 \t) — (16X + 2x i X> t ) V 3 (t) 

60 45 5 5 

4 (5) (4) / 16 \ (3) 
-—X' 4 V 33 {t) +12X 1 ' 3 T 3 (t) + [6U 2 + y$x - 12$ 2 J I 33 (t) 

(3) 

-(2f/X' ij ' - 2S(X' iJ ) + 4X{ iJ - 2Xf - 8P 1 , 3 ) T j (t) 

4 (6) 4 (4) 

+— r 2 x fc e mfcj ^'(t) + -X' k e mkj J m3 {t) 

1 (5) 9 (5) 2 (5) 1 (5) 

r 2 M^' fcfc (t) M iijk (t) +—x j M iijkk (t) M iij3kk (t) . (6.4d) 

15 W 15 W 15 W 30 W V ' 

The final term in the expression for S3 is purely a function of time, and as such does not 
affect the equations of motion through 3.5PN order. It comes in part from the surface terms 
Eqs. (CT), and in part from various integrations by parts of Poisson potentials to achieve 
the expressions shown. In A/3, all such terms cancel. Similarly, purely time-dependent terms 
which appear in A/3.5 and -B3.5 do not contribute to the equations of motion. 

As expected, the outer integrals make their first formal contribution to the field at 3PN 
order, however, the observable contribution vanishes to this order, so we have not shown any 
such contributions explicitly in Eqs. ( |6.4|) . In the next subsection, we study the contributions 
of the outer integrals in more detail, and show that through 3.5PN order, all contributions 
from the outer integrals are pure gauge terms. 



C. Outer integrals and the contributions of "tails" 

Our earlier qualitative discussion suggested that terms involving products of the 
monopole moment X and the quadrupole moment X 13 of the far- zone fields would contribute 
via the outer integrals at 3PN order. Because higher multipole moments involve higher pow- 
ers of 1/r or higher time derivatives, they would be expected to contribute at even higher 
PN order. Thus working through 3.5PN order, we might expect at most that products of 
X with quadrupole X K \ octupole X %3k or current quadrupole J %3 moments would contribute. 
Other terms, such as products of X with higher-order moments, or products of higher-order 
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moments, such as terms quadratic in X % \ will be 4PN order or higher. In studying the 
contribution of the outer integrals to the fields at 3.5PN order, therefore, it suffices to em- 
ploy the far- zone field given in Eq. ( |5.12| ). However, to illustrate the first non-trivial "tail" 
contribution, we will evaluate certain pieces of the outer integrals through 4PN order. 

We substitute Eqs. ( |5.12| ) into Eqs. fl4.4| ) using the "quick and dirty" rule expressed by 



Eq. (|4.9|) to determine which terms to keep, and obtain, in the far zone, 

A 00 = l4Xr- 2 n% ik (V h /r) - SZr^dyiF* /r) + KLr^tid^T 3 ' /r) - 2\Xr '// " T J 

•■14 8 

-llr-^d^ lr) - -lr- 2 h l d ljkl (V M /r) + -Tr~ l d l3k {T^ k /r) 

--Zr- 2 tid jk {T ]k /r) - 81r~ 3 h <lJ> d k (l ilk /r) - -Ir^n'd^V^ /r) 
3 3 

+ —lr- 2 n i e mb d ak {j bk /r) - 321r- 3 n <lj> e aki d k {j aj /r) - -Xr- 2 h i e akj d ik {j aj /r) 
3 3 

+0(pe 4 ) , (6.5a) 
A 0i = 81r~ 2 n j d ik (± jk /r) - %lr~ 2 Wd ]k { T k /r) - 6Ir- 2 n l d jk ( ± jk /r) 

( 3 ) 8 (3) 16 

+81r- l d j (Ti /t) - -Zr~ x d 3k {Z ijk /r) + — Ir' 1 e iab d ak {J bk / r) 

o o 

(3) (3) 

-81r- 2 n j (l ij It) + 21r- 2 h l {l" /r) + 0(pe 7/2 ) , (6.5b) 
A* = 2Xr- 2 h i d ijk {X jk /r) - Slr^Wd^T 3 /r) + 21r- 2 h i d l (l jj /r) 

(4) 2 8 
-Sir' 2 T li --l r - 2 h% jkl (l jkl /r) + -lr~ 2 h l d jk (l^ k fr) 

2 •■ 8 ^ 16 

+-lr- 2 n i d ij (l jkk /r) - -Xr^d^T" /r) - —iT^n^d^J^/r) 

o _ on (3) 

+ -lr- 2 h l e ak] d ik {J a] fr) Xr' l e aki d k {J ai fr) + 0(pe 4 ) , (6.5c) 

3 3 

(4) 

A*i = -8Xr" 2 T j +0(pe 3 ) . (6.5d) 

All moments X lJ , 2 t: > k , and J 13 in these expression are functions of retarded time t — r. 
Notice that the term kept in A tJ is actually of 0(pe 3 ) (4PN order) according to our scheme, 
however because it has 1/r 2 dependence, it will yield a 4PN tail contribution of a form which 
we wish to keep. 

We expand the derivatives and evaluate the coefficients £% L and £| L [Eqs. ( |2.34p and 
( j2.36| )l for each term, throwing away all ^-dependent terms. Terms with 1/r 2 fall-off yield 
integrals over Legendre functions Ql, as in Eq. ( |2.35p . The result, through 3.5PN order 
(and keeping all formally 4PN terms involving integrals over Ql), is 

( i-oo ( 4 ) 8 r oo ( 4 ) 

(N 3 ) C -m = Z\-m <ij> jf T j (t - rC)Q 2 (CR - 3 J i 1 jj it - rC)Q (CR 

+ _(n<y> - 2S 13 - 2S lJ lnr) V 3 (t) \ , (6.6a) 



3 

( S /-oo ( 5 ) S /-oo ( 5 ) 

(N 3 . 5 ) c -M = l[--n <t3k> J i T jk (t - rC)Q 3 (C)dC - -W ^ T' JJ (t - rC)Qi(CK 
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o (5) o (5) 1 

__ r ( 3 ^<«> _ 2S ij ) V j (t) + —(5n <ijk> + 18n^ ifc ) T jk (t) | , (6.6b) 

f ,00 ( 4 ) (4) 1 

(^3.s)c-^ = j i T 3 (t - rC)Qi(CR + 4^ X* j (t) J , (6.6c) 

r S /-oo ( 5 ) 8 /-oo ( 5 ) 

(K 4 )^ = j|--n<^>y i (t - rC)Q 2 (CK - g jf Z*' (t - rC)Qo(CK 

16 /- 00 (4) 

+— n <afe> e^' / jF fc (t - rC)Q 2 (CR 

3 Jl 

Ifi /-oo (4) i 

+y ^ jf J> k (t - r()Q (C)dCj , (6.6d) 

f /-oo (4) (4) i 

(£ 3 ) c _^ = x{-8 jf X" (t - OQo(CR + 8(1 - lnr) X" (t)| , (6.6e) 

r Q /-oo ( 5 ) qo r 00 ^ 

(5 3 . 5 )c-^ = l[ + -n l jf 2^' (f - rC)Qi(CR + jf J* (* - <)Qi(CR 

(5) A (5) 16 (4) 1 

+4r X™ (t) - -n l V 33 it) e ai3 n l J aJ (t) , (6.6f ) 

3 3 J 

/oo (4) 
T j (t - r()Q (()d( . (6.6g) 

Using the recursion relations satisfied by Legendre functions, we can establish the general 
formulae: 

fx {t -r()Q L «)d( = I(i ^ T) X(t-r) 

1 roo 

~2lTi I x ' (t ~ r OKkn(0 - <fc-i(C)K , 

X(t - rC)Qo(CR = X(t - r) - jf X'(t - rC)(Qi(C) + q (OR 

+ / X(t-r- s)ln(s/2r)ds, (6.7) 

JO 

where prime denotes d/d(, s = r(( — 1), X represents one of the multipole moments of 
the system (X u and higher), and we assume that, in the distant past the system becomes 
sufficiently "stationary" that as s — > oo, X(t — r — s) Ins — > 0. Since for a binary system 
that becomes unbound (r — > t> qs) in the infinite past (because of gravitational-radiation anti- 
damping, looking backwards), X in the worst case is proportional to (d j dt) A X %3 ~ mv 4 /r 2 — > 
mfg/s 2 , then this boundary condition is satisfied (see []6D] for a detailed discussion of the 
past behavior of binary systems whose evolution includes gravitational radiation reaction). 
Repeated use of these identities allows us to convert many of the integrals in Eqs. ( p.6|) 
into integrals of higher time-derivatives of the expressions, which are thus of higher PN 
order, plus residual terms that cancel many of the non- integral terms in Eqs. ( |6.6| ). It 
is also useful to expand the retarded time t — r — s about t — s, and to separate the lnr 
terms from the ln(s/2) terms in the integrals, leaving only terms proportional to X^ n \t) and 
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N c - 


-Af 


-A-f- 


K- 


-AT 


= O(e^) , 


Bl 


-AT 


= 0(e 4 ) , 


B c - 


-AT 


roo 

= -81 / 



/ °°lW(t - s) ln(s/2)ds. In the end, the only terms that remain at 3PN and 3.5PN order 
are 

( 4 ) O f oc (5) >, 

I u (*)--/ X" (t - s) ln(s/2)ds + 0(e 5 ) , (6.8a) 
3 jo J 

(6.8b) 
(6.8c) 

/•oo (5) 

8X / 2* (t - s) ln(s/2)ds + 0(e 5 ) . (6.8d) 
jo 

As these are purely functions of time, they do not contribute to the equations of motion 
through 3.5PN order. Alternatively, one can show that the terms in Eqs. ( |6.6| ) turn out to 
be purely gauge terms through 3.5PN order. In fact, by making the gauge transformation 
hy- v —>■ h^ u — — + ri^ v ^ a a (the linear transformation suffices to this order), with 

(A ( 3 ) q roo ( 4 ) o r°° ^ 

£° = 1 - T l (*)+-/ l u (t - s) ln(s/2)ds - -x ij / (t - s) \n{s/2)ds 
13 3 jo 3 Jo 

9 /-oo ( 6 ) 4 /-OO ( 6 ) 

+ -r 2 / X" (t - s) ln(s/2)ds + — / X™ (t - s) \n(s/2)ds 
3 Jo 45 Jo 

+-x i e ikj / ^ (t - s) \n(s/2)ds , (6.9) 
9 Jo J 

{/-co (5) 4 poo (5) 

-4x J / r j (t - s) ln(s/2)ds + -x* / X Jj (i - s) ki{s/2)ds 
jo 3 Jo 

+ — / V" (t - s) In s/2 )ds - -e lkl / J jk {t - s) In s/2 ds , (6.10) 
45 Jo 9 Jo J 

we can convert the outer integral contributions to h a/3 in Eq. fl6.6p to a form consisting of 
nothing but a 4PN tail term: 

Ifi r oo ( 7 ) 

(JV + 5) c _ A r = -— Xx ij / X <ij> (t - s) ln(s/2)ds + 0(e 5 ) , 
5 Jo 

= 0(6^) , 

(B i VAT = 0(e 4 ). (6-11) 

Note that, to this order, N + B = 2g 00 , and only the gradient of the term in Eq. (|6.11|) 
contributes to the acceleration, hence this term can be thought of as a 4PN tail modification 
of the Newtonian gravitational potential, or as a 1.5PN modification due to tails of the 2.5PN 
radiation-reaction potentials. This result is in complete agreement with the near-zone tail 
contribution derived by Blanchet and Damour |T5| using matched asymptotic expansions 
within the post-Minkowskian formalism. 



VII. DISCUSSION 

We have presented a method for direct integration of the relaxed Einstein equations in 
a post-Newtonian expansion, applicable to equations of motion and gravitational radiation 
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from isolated gravitating systems. As a foundation for future work, we presented a solution 
for the near-zone gravitational field through 3.5 post-Newtonian order in terms of Poisson 
potentials, together with a prescription for ensuring that no divergent or undefined integrals 
occur. In subsequent work, we will apply the near- zone results to the derivation of equations 
of motion for binary systems of compact objects through 2.5 PN order, and including 3.5 
PN radiation reaction terms. Work on the 3PN contributions to the equations of motion is 
in progress. 

The results presented here can also be applied to the gravitational radiation wave-form 
and energy flux from binary systems to as high as 3PN order beyond the quadrupole ap- 
proximation. It can also be used to discuss equations of motion and radiation damping of 
systems containing spinning bodies, as well as the structure and evolution of fluid bodies. 
These will be the subject of future work. 



ACKNOWLEDGMENTS 

We gratefully acknowledge useful discussions with Luc Blanchet. This research was 
supported in part by the National Science Foundation under Grant No. PHY 96-00049. 
The initial phase of this work was also supported by Fellowships to CW from the Fulbright 
Foundation and the J. S. Guggenheim Foundation. CW is grateful to the Observatoire de 
Paris, Meudon and the Centre National de la Recherche Scientifique in France, and to the 
Hebrew University of Jerusalem, for hospitality and support during a sabbatical year where 
this work began in earnest. 

APPENDIX A: STF TENSORS AND THEIR PROPERTIES 

Throughout this series of papers, we shall make frequent use of the properties of symmet- 
ric, trace-free (STF) products of unit vectors. The general formula for such STF products 
is 



[1/2] 

" ' (Al) 



n <L> = j:yr i2 \ 2l l _*f [n L ~ 2P 5 P + sym(q) 



where [1/2] denotes the integer just less than or equal to 1/2, the capitalized superscripts 
denote the dimensionality, / — 2p or p, of products of n l or 5^ respectively, and "sym(q)" 
denotes all distinct terms arising from permutations of indices, where q = l\/[(2 p p\(l — 2p)\] 
is the total number of such terms (see |T4] , |2"6"f for compendia of formulae). For convenience, 



we display the first several examples explicitly 

(A2a) 

i i 5 jk + n j 6 ik + h k 5 ij ) , (A2b) 

5 

hn ij S kl + sym(6)) + ^-(S ij S kl + S tk 5 jl + S u S jk ) , (A2c) 

// ' //■' - -(n ijk 5 kl + sym(10)) + — (n l S jk 5 lm + sym(15)) , (A2d) 

9 63 



h <ij> 


= n ij - - 








= n ijk - 




= n ijkl - 


<ijklm> 
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<ijklmn> fiijklmn 



1 ™, 1 



n J = n 



(^ife/^mn + sym ( 15 )) + (fiv«$w«S™» + sym(45)) 
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1 ( 5 y 5 ^ 5 mn + sym (i 5 )) . (A2e) 



There is a close connection between these STF tensors and spherical harmonics. For example, 
it is straightforward to show that, for any unit vector N, the contraction of N L with h <L> 
is given by 

^ Lfl<L> = (2[3l)H P ^ ' fi) ' (A3) 
where Pi is a Legendre polynomial. 



APPENDIX B: CANCELLATION OF ^-DEPENDENCE BETWEEN INNER 

AND OUTER INTEGRALS 



Here we demonstrate explicitly the cancellation of ^-dependent terms between the inner 
and outer integrals. We assume that, at each iteration step, from just inside the boundary 
of the near zone out into the far zone, the source stress-energy tensor 7v-iA a/3 can be decom- 
posed into terms of the form fB,L{u)n <L> r~ B , where u = t — r is retarded time, and n <L> 
is a STF product of unit radial vectors. We calculate the behavior of the inner integral of 
such a term as the integration variable approaches TZ from below with the result obtained 
from the outer integral of the same term. We consider far-zone and near-zone field points 
separately. 



1. Far-zone field points 

The inner integral is given by Eq. (|2.13|) , with the multipole moment given by Eq. ( p,14|) . 
We want to examine the behavior of the moment, as Ix'l —*Tl, that is 



lo7r J r 



l<L> 



7 E ^/SWG^q^A^ , (Bl) 

4 m=0 m - 



where the superscript (m) denotes m retarded time derivatives, and where 

A L ><2 = — I h <L> h Q dQ , (B2a) 
An J 

im (T)\ I „/ 2 +9— B+m 



Then, from inside 1Z, 



ft 3+9 "- B+m /(3 + q- B + m) 3 + q - B + m ^ 
In ft 3 + q- B + m = 



(B2b) 
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oo / -i\q oo / -i \m / 1 \ 

h^B,L - E ^# E [ ~^rd Q ( ) G™ Q (K)A L >$ . (B3) 



It is straightforward to show that the contraction of 8q with A L '® is given by 

q < L 

A L > Q dn = { L + q = odd . (B4) 

2^l((g+L)/2)l | V 2|(q-L)/2g = > £ 

(g+L+l)!((g-L)/2)! I V I a <^> « ^ ^ 



Using the fact that 

T*(fM\ f 



r 



(B5a) 



(see eg. [Q) and redefining summation variables, q = m + q — k, j = L — k, we obtain Eqs. 
( [2.371 ) an d (P-38|) . Evaluating the outer integral for the same term yields ^-dependent or 
In ^-dependent terms that are precisely equal and opposite those of Eq. (|2.38|) . 



2. Near-zone field point 



In the near zone, for |x'| > |x|, Eq. ( 2.22Q together with the specific decomposition of 
A a/3 gives 



h°0 



-r oo 

— E 



E 



-X 



m=0 



ml 



n 



Q 9t fBAt 



n 



l<L> 



& Q {r' m ~ 1 )r' 2 dr'dQ! . (B6) 



We use the fact that |JH | 



ml 



^ (2q-4k + 1)!! 
^ (2g-2A; + l)!!(m-2fc)! 



(m-2Jfe-l)!! 
[m-2q + 2k- 1)!! 



x 



2 k k\(q-2k)\ 



§K^I<Q-2K> r im^q-l 



(B7) 



where k m = lesser of {[g/2], [ra/2]}, 5 K denotes a product of K Kronecker deltas, the quantity 
in square brackets can be evaluated for negative or positive values of the arguments, and the 
expression fi K f l '<Q- 2K > j s to be symmetrized on all indices (since the expression ultimately 
is to be contracted on x® no explicit symmetrization is needed). It can then be shown that 



f,Q— £> 5 K fi' <Q - 2K> h' <L> dVL' = — 5r - oifi 

An J (2L + l)\\ L ' q - 2h 



<L> 



(B8) 



S^ =0 (— l) n (t)r /n /n\, integrate over r' toward 1Z, rearrange 
the summations, and define r = (q — L)/2, and q = m + n, and obtain 



We then expand fit — r') — ^ n=0 
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1 tt/J 
N^M B,L 



ry \ B—2 OO 

ry 



5=0 



with 



[-l) L+ i{2) 2+L - B {L + r)\ 



[9/2] g 

BM ' f^ m t^ r (q - m)\(2L + 2r + l)!(m - 2r)!r! 

z q-L-2r-B+2^ q _ L _ 2r - B + 2) 

In ft 



X 



rf g /B,L(t) 

dti 



(B9) 



(m- 1 -2r)!! 
(m - 1 - 2r - 2L)H 

q-L-2r -B+2^0 



q — L — 2r — B 



0. 



(BIO) 



Here too, evaluating the outer integral for the same term, for each B, L and q yields in- 
dependent or In ^-dependent terms that are precisely equal and opposite those of Eq. ( |B10| ) . 



3. Source terms with lnr dependence 

Until now we have assumed that the stress-energy source A af3 can be decomposed into 
terms of the form fB,L{u)h <L> r~ B . At sufficiently high PN order, tail contributions to the 
fields will arise, leading to the possibility of lnr dependence in A Q/3 . To illustrate that 
cancellation of 1Z dependence occurs in this event also, we consider source terms of the form 
f B ,L( u ') fl><L>r '~ B ' [nr '- Noting that, from Eq. ( gggfl , lnr' = -ln[2(( - y)/r(( 2 - 1)], and 
incorporating this logarithmic term into the outer integral, Eq. ( p. 27] ), we obtain 




2 J-i ^' "Jay) V(( 2 -l)J V(( 2 -l), 

xf B , L (u-r{C-l))- 
d 

= ~dB Nk C-^B,L- ( B11 ) 

For the inner integral, the only difference which the logarithmic term makes is in the radial 
integral, now given by 

r'R 



G b,l,q W (lD) = / r' 2+q ~ B+m In r'dr' 
d 



dB 



Gl L ATZ) . (B12) 



Thus, if the original coefficients cancel for all B (and if we can treat B formally as a 
continuous parameter), then the coefficients generated by lnr terms cancel. 

An alternative method is to show directly from the definitions (eg. Eqs. ( [Bll|) and 
( |B12[ )) that, for both the inner and outer integrals and for z < 1 and z > 1, 



N h a B P L ln) = In ft N h a B p L - J N h B p L dz/ z , (B13) 

modulo z- or ^-independent terms. Then, if the ^-dependent parts of tv^bl cancel between 
outer and inner integrals, so too do the ^-dependent parts of n^bl- 
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APPENDIX C: BOUNDARY TERMS 



The boundary terms in h°$ that arise from integrating by parts various integrals over 
A4 are given by 

N dM =A<f> T°i(t,x')d 2 S' + lr 2 dU ^(t^dPS'.-U'dtl T ij (t,x.')d 2 S' 

JdM o JdM o JdM ' 

--x'd 2 I (t, x')x ,l d 2 S' + —r 4 df I (t, x')d 2 S' 
3 * JdM 30 * JdM K ' 

_A r Va 4 / T°3(t,x')x' i d 2 S'-—r 2 x i d?<{ r ij (t,x')d 2 S', (Cla) 

15 JdM 15 JdM 

Kl M =Al r ij (t, x!)d 2 S' + lr 2 d 2 t I r ij (t, x!)d 2 S' + \x k d\ I r 0j (t, x>)x nh d 2 S' 

JdM J 3 JdM J 3 JdM J 

--x k d 2 I r^(t, x')x ,k] d 2 S' - -df I r°i(t, x / )r' 2 x' i d 2 S' , (Clb) 

3 JdM 9 JdM 

B% M = -Ad t £ r k ^(t, x!)x' ]) d 2 S' k - 2d 2 I r ok (t, x!)x' %3 d 2 S' k 

JdM " JdM 

JL r ^dl I r k ^(t, x!)x' j) d 2 S' k - -r 2 d 4 t I r ok (t, x!)x' ij d 2 S' k 

3 JdM 3 JdM 

+ -x l df I r^(t, X ')x" l) d 2 S' k + -x l d? I r ok (t, x')x' ijl d 2 S' k 
3 * JdM 9 * JdM 

+ lx l df I (r"*(t, *!)x' l]j + r k ^(t, x!)x' l]i )d 2 S' k 

9 JdM 

/ rjk/j x ')( r /2 x /ii) x + f ok, t x 'yV l V<% 

18 JdM 

+ -d? i (r k V(t,x')x ,l]jl + r k \t^)x' m )d 2 S' k 

3 JdM 

-—r 4 d" I r k( *(t, ^)x' ]) d 2 S' k - — r 4 a 6 / r ok (t, x!)x' ij d 2 S' 
30 JdM 60 JdM 

+—r 2 x l d b t I r^Ht, x')x' jl) d 2 S' k + —r 2 x l d\ I r ok (t, x')x' ijl d 2 S' k 
15 JdM 45 JdM 

+-irVdf / (r^(t, x.')x ,l]j + r^'(t, x!)x' l]i )d 2 S' k . (Clc) 

45 JdM 



APPENDIX D: PROPERTIES OF POISSON POTENTIALS 



Here we list some useful properties of Poisson potentials and superpotentials, given by 
Eqs. (|4.10j ). These rely upon the general result, which can be obtained by integration by 
parts, 

P(V 2 g) = -g + B P (g), (Dl) 
where Bp(g) denotes the boundary term, given by 

n 2 dtt' . (D2) 

r'=Tl 



Bp(9) = -r i P^9;in(^(t,x')|x-x'|) 

47T JdM X — X 
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The boundary terms must be carefully evaluated case by case to determine if any 1Z- 
independent terms survive. All 1Z- dependent terms can be discarded. Some useful formulae 
that result from this include: 

P(\Vg\ 2 ) = -\{g 2 + 2P(gV 2 g) - B P (g 2 )} , (D3a) 

^(Vs ■ V/) = ~{fg + P(fV 2 g) + P(gV 2 f) - B P (fg)} , (D3b) 

P(f\ W| 2 ) = ~\{fU 2 + P{U 2 V 2 f) - 2E(fU) + AP{UVU • V/) - B P (fU 2 )} . (D3c) 
In many specific cases, the boundary terms can be dropped: 

(D4a) 
(D4b) 
(D4c) 

, i2 _ '*»-, (D4d) 

' ij ~ Y2 Y ' ij + \ xky,l3k ' ( D4e ) 

while in others, there are contributions from the boundary terms. For example, in the 2PN 
potential P(VU ■ VX), the boundary term yields the term | J M a(t, x)d 3 xd 2 J M a(t, y)d 3 y. 
Using Eq. (|5.6p, we obtain, to the necessary order, 



P{U) = 




P(X) = 




P(\VU\ 2 ) = 




P(x i U^ k -) = 


--x l X 
2 


P(r 2 U' 11 ) = 


--r 2 X 
2 



1 



1, •• 1 



P{VU ■ VX) = --{UX - E(X) + 2G 2 - -ZT\t)} + 0(e 5 ) . (D5) 
Similarly, we find for the 3PN potential, 

(4) 1 (4) (4) (4) (4) 

P(W -VY) = --{UY + 12P(U X) - 21 T\t)} + 0(e 5 ) . (D6) 

For the Poisson Superpotential S(f), we have 

S(V 2 g) = 2P(g) + B s (g) , (D7) 

where 

B s {g) = ^i U, xOlx-xKlnff^^l K W . (D8) 

47T JdMl \ |X — X'| J \r'=Tl 

Thus, for example, in the superpotential (d/dt) 2 J M r 00 |x — x.'\d 3 x', we find the term 
S(V 2 U 2 ) = 2P(U 2 ) - 3(d/dt) 2 ( [ ad 3 xf + 0(e 5 ) 

JM 

(4) 

= 4Gi + 4G 2 - 31 T\t) +0(e 5 ) . (D9) 

Other useful identities include 

= x'U-X'*, (DlOa) 

E(x ij ) = -Y' ij - 5 lj X + x lj U - 2x {l X' j) . (DlOb) 
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